{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Bowers method with well log\n",
    "\n",
    "Pore pressure prediction with Bowers' method using well log data."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Predicton of geopressure using Bowers' model needs the following steps:\n",
    "\n",
    "1. determine Bowers loading equation coefficients A and B\n",
    "\n",
    "2. determine Bowers unloading equation coefficients $V_{max}$ and U\n",
    "\n",
    "3. Pressure Prediction"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {"nbsphinx": "hidden"},
   "outputs": [],
   "source": [
    "import warnings\n",
    "warnings.filterwarnings(action='ignore')\n",
    "# for python 2 and 3 compatibility\n",
    "# from builtins import str\n",
    "# try:\n",
    "#     from pathlib import Path\n",
    "# except:\n",
    "#     from pathlib2 import Path\n",
    "#--------------------------------------------\n",
    "import sys\n",
    "ppath = \"../..\"\n",
    "\n",
    "if ppath not in sys.path:\n",
    "    sys.path.append(ppath)\n",
    "#--------------------------------------------"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "from __future__ import print_function, division, unicode_literals\n",
    "%matplotlib inline\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "plt.style.use(['seaborn-paper', 'seaborn-whitegrid'])\n",
    "plt.rcParams['font.sans-serif']=['SimHei']\n",
    "plt.rcParams['axes.unicode_minus']=False\n",
    "\n",
    "import numpy as np\n",
    "\n",
    "import pygeopressure as ppp"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1. determine Bowers loading equation coefficients A and B"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Create survey with the example survey `CUG`:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "# set to the directory on your computer\n",
    "SURVEY_FOLDER = \"M:/CUG_depth\"\n",
    "\n",
    "survey = ppp.Survey(SURVEY_FOLDER)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Retrieve well `CUG1`:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "well_cug1 = survey.wells['CUG1']"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Get velocity log:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "vel_log = well_cug1.get_log(\"Velocity\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Preprocessing velocity data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "vel_log_filter = ppp.upscale_log(vel_log, freq=20)\n",
    "\n",
    "vel_log_filter_smooth = ppp.smooth_log(vel_log_filter, window=1501)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "View velocity and processed velocity"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAARMAAAHyCAYAAADIsSqfAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAIABJREFUeJzsnXd8VFX6/993SjpJID2BACHkhq5EihQBkaLCWsC6gqJiXfu6q+sq637VXV3dVdfys6G7ClhWERBFEaWDQAApSS4YaiAJCZCeSabc3x+TuWYyk2SSzGQSct685pW5957y3AzzybnnPOd5JFVVEQgEgrai87cBAoHg3ECIiUAg8ApCTAQCgVcQYiIQCLyCEBOBQOAVhJgIBAKvIMSkAyPL8gRZlv/tw/ajZFn+uQXl9/rKFg/7b9XvQ5blJ2RZ/lGW5Y2yLPeqO/c7WZZ3yLL8lvct7ZoIMen4+MQRSJblEOADIMLftrSQFtkgy3IcEK0oyiTgn8C9siz3BH4LjASOyrI8xftmdj0M/jagKyLL8gJgvaIoP8qy/Db2L/VTQDjwtqIoH7SgrSuAx7D/YbgHKAZeUxRlpizLE4DZiqLcJ8vyp0Bv4CiQCbwCzAXWNdP+COBN4AgQUnduFPASIAF/BpKBQGAssByYDzwEvAacAi4A5imKslGW5cVAInAGmF1X70MgAdijKMrdsiyn1J3TASsURXmuEdvGAs877l1RlN0NyyiKUlhnC0AMUIhdRFYpimKTZXkNcCmw2k37QcDnQBiQoyjKnbIs5wMKUAQMBO5VFGVtU7/DroIYmfiHpcDUuvd9gGuBd4BxwP2eNiLLsgT8FZgA3ID9SwLOf71VWZYjsX/hxwDnK4ryvKIoJkVRznrQzV/r7HsA6FF37l/A1di/hC8BOwEZu3CkA3qgqs6eh7AL5VWyLPcA+iqKMhG70HQD7gDWKooyFjgsy3IMdrG5G/gNcFcTtr2CXZDmAguauglZlsOAecB72EXxZN2lM0BsI9UGAuWKokwAvqn7fRuAPwCRwL+BAU3125UQYuIHFEXZA6TJstwf+xcxDfuXbg0QWPfl94QYoEBRlFpFUQ4BDZ//A+t+VmD/cm/CLgQtIVFRlEOKopzA/tcYIFRRlFOKopQBJuyjljSgBuiJXVQAtiqKUoD9i2tUFOUM8Lksy6uAWXU2ydhHSgD/UBTF0cfz2IWqqceavsBi4G0guJn7eAV4ts7minrlQ2n8e7AL+2PQKuACRVFUoACwAHmADfvoTIAQE39yCJgDfIF92PxI3XP9y9i/lJ5QBMTLshwgy3IfYAVg5tcRxMS6nwOBzYqijFYU5Y0GbTT3ZSiQZbmvLMsJQHTduXJZluNkWY4AAhRFKccubEewf8Edk7qm+g3JspwE5CmKMh0IwD4SywEy6oosrHus+iP2v/6/p+lH8Z+By+peSxsrJMvyfcAZRVFW1J3KxD5KAxgMHGuk6lBgdZ29o+oevwSNIMTEf3wJ3KAoyjbgBWCBLMsbgQGKolR70kDdX8oF2Oc9/gM8oSjKSaBUluXX+XVy9QhwkyzLa2VZ/lSW5UH1mmluQvNp4H/Y500cj0UPYxfBr4FH6s7tAX6p68sxd9Gw7UJgtizLm4B+wA7gXWCSLMs/AhWKomzHPu+yCFgI1NaJljsWAN/W3X+5uwKyLIcDL2IXgx9kWV6gKMpxwCDL8mvAE8DHjbR/FHhUluXNde0frXdPHWEyukMhiV3DnQNZlqdhn2hVsY8mVOzD9jUe1B2Pfd7CWvf6R8NJw7pVj48btL9SUZSXvHgbPkWW5ZuBW3C+h/sURdnvpqwemAHsUxQlV5blJUBc3WUJ+whqTrsYfo4gxEQgEHgFvy8Ny7Ici32IXgbMr5sgEwgEnQy/iwnwaN0rEPsS32v+NUcgELSGjjABG68oyj5gL/aVAIFA0AnpCCMTx/KhHvuav1syMzPF5I5A0A5kZGS0ynemI4jJPlmWM7A/5hxvqmBGRkZTl9tEdnY2Awb4zpmxs7ffHn2I9v3fR2ZmZvOFGqEjiMn7wOvY9z/c4WdbBAJBK/G7mNSt3oj1fIGgk9MRJmAFAsE5gBATgUDgFYSYCAQCryDERCDwIfPmzQOgvLyc++93DVWzdOlS1qxpdnuVxtatW9myZYvX7PMmQkwEAh8yYMAAcnNzyczMZPTo0W1ub/To0Vx44YXa8eOPP05FRUWb2/UGQkwEAh9y4YUXsn37dnbs2EFaWhp/+tOf+N3vfseJEydcyubk5HD//ffz0EMPUVhYSFZWFg888AD33HMPhw8fBuCHH37QRjJLly5l3759vPzyy+Tk5PD++++za9cuCgsLee45t5EufYrfl4YFgo5ATk4OBQUFra4fHx/v9nxGRgZff/01lZWVLF68mISEBAwGA9nZ2SQlJTmVfeedd3jhhRfIz8/nww8/pKCggOeff57Kykq3o4+rrrqKbdu28eCDDxIWFkb37t1ZuHAhPXv2ZMaMGa2+l9YixEQgANLT00lPT29TG9nZ2S7nQkJCqKysJDw8nKqqKu655x6KioooLS11KStJEjabDVVVnV4Wi4WDBw/Su3dvlzp6vR6r1QpAXFwcZ86coaysjDlz2t91S4iJQOBj0tLSSEpKIjU1lcceewybzcZTTz3lUu7WW2/lscceQ5IkHnvsMc6cOcMf//hHzGYzv//97922fdFFF/HEE08wadIkZs2axciRI9s0wmoTDVWwo7527Nih+pKsrCzRvp/7EO23rY9ffvlFveuuu9SCgoJWt1/3PWvVd1SMTASCc4R+/frx5ptv+q1/sZojEAi8ghATgUDgFYSYCAQCryDERCDwMY8//jgPPvggd955J7m5uX6z49577/Vp+0JMBIJ24JlnnuGpp57ik08+8ZsNkuTbTKZiNUcgaAdUVaW8vJw9e/bw2GOPYTKZeOqppzh16hRvvPEGer2exx57jNOnT/PWW29hNpt59NFHNc/ZsrIyHn/8cU6cOMGrr75KQEAA9913H1FRUfzpT38iJCSE2bNnM2jQIF544QWqq6uZMWMGgwcP5sknnyQhIcHn9yjERCAABr8xmP1FLon/PGZQzCA+m/RZo9effvppDAYDc+fO5euvv+a11+wZXZ599tlGXejLy8tZuHChkwv+mTNnCAsL49prryUxMZHi4mJsNhuzZs1i4MCBrFixApPJRHx8PPv27eP48ePcfPPNpKWlsWDBglbfnycIMREIgH337GtzG+7c6R385S9/ISwsjG3btjFy5EjtfFMu9L/88guAkwt+YGAgt956K8uXL+f48eNMmjSJBx98kC1btrBjxw6io6OZOXMmY8aMYfXq1RQXF2MwGDAajag+zt4pxEQg8CPNudDPmzfPyQU/Ly+PZcuWYTKZmDx5MrW1tbz77ruoqsqUKVMYO3YsTz75JJ9//jnXXnstGRkZPP/888TExPh8zsTvbvKevoQ7vX/bb48+RPv+76Mt7vRiNUcgEHgFISYCgcArCDERCAReQYiJQCDwCkJMBAKBVxBiIhAIvIIQE4FA4BWEmAgEAq8gxEQgEHiFLudOb7FYMJvNmmux46fZbMZkMmnl1Hr7GFRV5VjpMYx6I/Gh8Y2Wcfe+W7duvndjFgg6AF1KTGpra/nb3/7Wqrqb2Uw++cxiVovrPvroo63qUyDoTHQpMQkICGh0G3ZWVhaffdb4FvLhDOdd3mUNa7iYi5H4dbQxZ84cUlJSvG6vQNCZ6FJi0hSSJLkVmuXLl3P69GmOHTvGPObxIR9SSy23xN+ChITFYuG7777j9ttvx2AQv05B10VMwDbDyJEjmTBhAjfccAPJUcn8IfYPlHYrZWfwTvr3709gYCCFhYU8++yzPP300z6PGSEQdFTEn9JmqJ+QOi0tjYMHD3Js8THeKn+LoMNBRBEFQEREBOeffz7btm1DlmUiIyP9ZbJA4BeEmLSQ/v37c8OlN3BkwxHWm9dzT4970Ov15OXlsW3bNqqqqggLCxNiIuhyCDFphpqaGvbu3UtISAhmsxmz2cyZM2e4OvFqHjn8CEelo8hGGYCqqiouuugiBg0a5GerBYL2R4hJIzz99NNOx+PGjSMoKIjAwEB69uxJamoqCxIWsDBrIfdefi96vV57lZaWau91Oh16vV7MpQjOeYSYNMIf/vAHTCYTVqsVq9WKxWJh3759fP/991qZSirJkXJYsWIFNptNK2u1Wp2OHe9jY2O56667hBOb4JxEiEkjbNy4kZ9++glJkpxewcHBSJJEVVUVwQRjU20kDU2iX2Q/9Hq9JjxmsxmDwcCQIUOQJIns7GwGDBjg79sSCHyGEJNGmDJlClOmTGn0usVioaKiAsNWA/dtuY/n+j/HgZ0HXMotXboUgJSUFCEmgnMaISYtRFVV7fHFZDJh+8lGb3rz8K6HmctcQghxKj958mTGjRvXZE4VgeBcQIiJGyoqKnjppZeaLNO3b1+MRiODBw9mmH4Y7+a/y07bTm6Lug29Xo/RaESSJNLT06mpqdHmTnQ64ScoODcRYuKGsLAwt671VVVVnDx5kpCQEJcJ1r/0+wvTVk5jeNVwIi2R1NbWArB7924AdDodn3/+udbWddddR3p6evvckEDQDggxaQEhISGkpqa6vXbs2DHOrzmfz/iMm7gJPXrt2sSJE4mNjRVzJoJzGiEmbaCoqIg33ngDsO9IHsMY8g357Evax1+G/4XU1FRCQuxzKGLORHCuI8SkDcTExPDUU09hNpspKyvj9ddf5wrLFXxw9AMeOPoAE5nott6UKVMYM2ZM+xorEPgYISZtRJIkSktLKS8vJyAgAGrhRm7kbd4mmWRSSEGSJMLDw7n33nsxGo3+Nlkg8AliacEL7N+/n40bN5KUlETfvn0Z2nco82Pns1K/EhMmVFWltLSUwsJCf5sqEPgMISZeYOLEicydO5e5c+cSGBjI4cOHCT8VTh9rH1axSiv33nvvcezYMT9aKhD4DiEmXqS0tJSjR48SFBQEwFSmcoQjHOAAkiQxYcIEkpOT/WylQOAbxJyJF3HMi0iSxMqVK8nNzWV6zXR+4Af6q/1Zt24d69atIyUlhTlz5vjbXIHAqwgx8SKSJBEaGsovv/xCVlYWAGmksYY15JLLnZPvZOjQoeh0OiorKwkICBATsoJzBiEmPiA/P197r0PHOMaxkY2krkll27Zt2Gw2bDYbffr04dprr/WjpQKB9xBi4mVKS0sZPHgwgwYN4tixYyxbtoyBDORLvkRFZdy4ccTFxdG7d29/myoQeBUhJl5EVVUtNYYj/glACSVEEomExDfffKOVNxgMPP7442Lzn+CcQIiJF5EkyWViddWqVSz8aSGxxLqUHzFiBNXV1RgMBiRJ0kI9CgSdESEmPsRms5GZmYkRIxYsLtf379/P7t27sVqtqKrK6NGjufjii/1gqUDQdoSYeIldu3axfPlyt9ciiaSEEpKTkwkMDHQKNq2qqvYqKiri448/5qqrriIwMLCd70AgaBtCTFqJ2WzmhRdewGJxHXE0xIYNCclj79e///3vgNgQKOhcCDFpJUajkSeeeMLttaqqKmpqarDZbJw6dYpHP32UJJJcyt15551OGQMFgs6MEJNWoKoqP/74I2azGcAlgn15eTk///yzVr6IImKIcWnn/fffd1r1qf/Io6oqkiRxxx13EB0d3T43JhC0ASEmreTw4cPk5eV5VFZGZgUruIALCCJIO+8I7eiOKVOmkJGRYQ9rIBB0AoSYtAJJkpg5cyb5+flOmfpUVcVisWCxWLBaraxZswaAfvQjlVSWs5xruAYJ90m4goKCkCQJg8HAiBEjhKu9oFMhxKSVxMbGEh0dzf/93/+h0+mw2WzataSkJAIDA5FlmdraWk6ePMkM6wzesr7F4R6HefySx53Siep0Os1xzfGIU15eTlhYmBiZCDoNQkzagE6n46mnnuLzzz9n//79REVF0a1bN0pLSzlx4gRRUVGoqkpISAg1Z2uYxSzeO/0ecZ/EEUqoS3vdu3d3Op42bRqyLLfX7QgEbUKISRuRJImIiAiioqKQZZno6Gj0ej2SJGGz2VBVFZvNRl5eHrt27WIQg/iO77iSK10ed86ePet0/PHHH3tsR1BQELIsC9d8gd8QYuIFHKlEKysrKSkpwWAwaA5pS5YsoaSk5NeyTGEhC9nKVi7kQo/a1+l0pKamoqoqUVFRTJs2zaVMdna2EBKBXxFi4kWOHz/O2rVrG431GhMTg9Fo5BEe4f8K/o/JPSfTK7iXNm/ScInYQXBwsJavx5E6QyDoaAgx8SLp6emkp6dTW1vL8uXLsdlsTvlyioqKtPdDGMKKYyuYzOQm2xw8eDChoaFi7kTQ4RFi4gMMBgMDBgzAYrFQWFjIsGHDOHnyJIqiaGXSSWc5y5nEJHSNhOINDg7myiuvFDuJBZ0CISY+wJGw3GazkZ6ejtVqJScnh0svvZSzZ8+yZMkSetb2JIQQNrGJ8YzX6l5++eVccMEFfrReIGgdQky8yI4dO1i5cqXTuYCAAHQ6HSaTyem8Dh2zmMXbvE0KKdrenZUrV7Jq1Spt/gTsPi2333670zmBoKMhxMSLXHDBBQwfPhyLxcLJkycJDQ3VVnV++uknevTogdVqxWq1snbtWiKIYCADOcpRp42ARqORyMhICgoKADh58iR//etfG+23d+/eJCcni8ToAr8ixMTLHDhwgM8//1wLTdCjRw8kScJsNmsxSuqv1BznOIMY5NSGyWTShMRBamoqwcHBgP0xavDgwRgMv3589ZefBQJ/IMTEy6Snp3PTTTexZ88epx3BeXl5LgJRSSVnOUsvejXb7i+//OJ0PHLkSCdxKiws5NSpU4B9CVo8EgnaGyEmXsZqtVJVVUVQUJC28c9qtboICUAVVYQRhp6Wr9a8++67TV6fPHky48aNa3G7AkFrEWLiJSorK3nxxRedzo0fP57ExER0Oh3BwcH07t3byUFtZ+FO1v20jluvvJVevZofnTRFdna2mDMR+BUhJl4iNDSUBQsW8P3337Np0yYANmzY4FTGcd5BQUIB3Szd2LJlC1u2bKGmpoZDhw5p1x3zLZIkaQJU//Fl8uTJmmesQOBvhJh4mUsuuYRLLrlEOy4tLeXll192W7Y8v5wDHGBb2Ta60c3l+pkzZ5rsa9GiRYB9FSklJaUNVgsEbUeIiY9YtmwZu3fvbrJMKqkMZSjv6d/j1qhbGRI8BICjR4963E9ISAhZWVns3buX/fv3M3v27DbZLRC0Fp+JiSzLc4BI4BPgP0AZMB8IaupYUZQyX9nkS1auXMmOHTvcXqsfMc1qtWrzJnq9nst1lzNQHcj7xe8zLGAYV3Wzp7moqalpts/k5GT69+8PwKlTp0hISHASsPpL0PVxnA8MDGTgwIFi5UfgFXwiJrIsDwUmAHuAR+tegcBcoHczx6/5wiZfc/nll3P55Zc7BYQ+efIkJ06cQFVVTp8+TW5uLiUlJU5R2QBiiGE+81lhWsFLppeYy1wMHnw0x44dc0qfsXfv3hbbHRkZyf333y8ERdBmfCImiqLskWX5v8D5QLyiKPtkWQ4AbgRimznu1NSfJDWZTHz77bce1QshhGu4hk/4hFWsYgYzWtTvtGnTGD16dIvtFQi8hS/nTBx/6hybUvRAWDPHrrEM61F/O7+3MZlMXmlfVVUOHjyopcEYOHAgABaLRdv9q6oqOTk5LnV16LiKq3ibt9nDHoYy1ON+v/32W024Zs6cqXnLehNv/Y5E+x27j9bSHhOw+2RZzsD+GHMMqGji+HhTDfnSj8JbfhqOYNCVlZVOjw7FxcVERUW5LBc3JIggruVa/st/iSfebcLz6Ohozj//fHr27AlAr169yMnJ8bmfia99WUT7/u8jMzOz1XV9LSYq8D7wBvZRyB1ADfB6E8edGkmSGDVqFAAfffQRubm5LW4jnnjGMY71rGc2rqszxcXFrF69utl2goKC+P3vfy/iofiZmsrdWEy/YLNVYghIIih8ElVnlhIUcTF6Qw9s1kqqS75B0gUSGDYGvTHK3ya3Cp+JiaIo64B1dYdzGlxu7vic4KabbgLgm2++Yd++fRgMBm1ytqqqymUitj7hhHMMz3ITN+TGG2/UVnkE/icw9Dz0higsNUcICh+PqXwjKr9+9jUVmwkKH4/OEEV16RpCIl1j/HYGRATidkCWZaqqqigrK6O8vJyKioomhQRAQSGF1jmi5eTksGbNGqxWa6vqC3xLULdx6PTh2rHNWoHeGIckGZCkzuv61Xkt70SkpKQwfvz4ZudLHJgwcYADTGd6q/rbuXMnYN9prKoqN9xwAxEREa1qS+B76guIajM1UbJjI8TEx3jiCduQ/eynL33dJupqjlmzZhEZGemUSD0sLKzF7QjaD50hBkvtCfSGaFQ672hSiImPueKKK7jiiiv44IMPKCoqQpIkKisrm6yzi11OcWEbMmOG3QclLCxMi1ovdg13XgJDz6O67AdUWy1B3TzLpdQR6VRiYirf7LO2bVajT9s/b5hMcmIfux9KbhHffu/eV6CIIqoDKnj7nrsw6Br7eGq0nw6bfW1/e/Rxzrcv6bXrOn045uoczHWXdPpI0IPZlIvZ1PgKoO8/58BW1+xUYhLUbYzP2tbpswnq5ru/7IFB2WzbdZSffvqpyXJHOEJSbTIvvfyj03lJkrRE5zfddBPh4eHasSRJSLpcgroN9Jn94PvfkWi/I/TRcf1MBPWYPn0606f/Oqm6Z88eli5d6lQmhRR+4AdOcIIkktDr9VitVi1qm8Vi4bPPPqOszHU/5F133UVcXJzP70MgcIcQEz+gqio7duzg66+/drkWRRRXcAWLWcyN3EiS9deo9Xq9nmuuuYZ+/fo5BUwC+5yJEBKBPxFi4ifcCYmDdNIBWMQi5jGPGGIAe/iCjz/+WCvXp08fbQ9OeXk5+/bt066FhoZy6aWXit3AgnZDiIkfkCSJBQsWOJ07c+YM33//PZIkkZWVRTrpFFLINrZxOZdr5fr37094eDiqqmKz2bBYLISFhREZGUliYqJWLjAwUAiJoF0RYtJBCAkJYfDgwVpK0S+++IJooskn36ncwYMH3dY3GAz88Y9/dMqlIxC0J+J/XgfBYDDw2WefAWgb8yxYMGJ0KdutWzfCw+3u2FOnTiU5OZns7GwhJAK/Iv73+ZmCggKOHz+OqqpMmTKF1atXa3tqyiknjF+9V6+77joXb9aYmJh2tVcgaAwhJn6muLi40cnY05x2ykHct29fLYufQNDREGLiZwYMGMDIkSPZtm2by7VqqgkhRDv++9//DsDQoUO56qqr2s1GgcAThJj4Gb1ez6WXXsqll16KqqqYzWYt+XkkkZzmtFZ22rRpjBo1SqzSCDokQkw6EHv37nXyiE0mmZ3s1Db9OeK8jhs3jsmTJ/vLTIHALUJMOhB9+vRxPqYPy1iGFatTcvONGzeyceNGp7LXXntte5goEDSKEJMORHh4OAsWLCAnJweTyYSqqry34j3OqGc0L9iGjBo1ij59+jSacEsgaC+EmHRA0tPTqaio4OWXX6af2s/JC7ah56yDjpr+QNB1EGLiB44ePcoHH3zgUdlJTOIN3mAgA+lLX9577z23E7BVVVWsW7cOi8XC7bffTlBQkJetFgiaRoiJj7FYLPzjH/+gtra2VfWDCeY3/IYv+ZL7uZ+8vLwmy995553CF0XgF4SY+BiDwcDjjz/eZFjFt99+m/z8fLfXAGKJxYQJCecRyYABA4iKitJyGcfGxhIbGyuWjgV+QYhJB+COO1xzj1VUVPDhhx9y6tQp9rCHIQxB1yAzScN5kpycHNavX++2jyeffBKdTmQ2EfgOISYdFIPBwKlTpwB7XNjW5tAxGo387ne/E0Ii8DlCTNoZVVWxWq1UVFTwyiuveFTHgsXJz6Q+Y8aMISAggKKiImJjY4mLiyMxMVF71JEkidDQlqfMEAhaihATH6KqKl988YVTBLSWUkophznMJVwCQL9+/YiPj0en0xEREcF5552HXq8XqS4EfkeIiY+Jiopi5MiRnDlzhqioKKfkWO5e1dXV7NmzRxtZLKpexGhG04MeAOTm5jolQ//qq68YOnQoPXv2pKSkpNn2649YwP44JR6BBN5AiIkPkSSJiRMnUllZSU5ODpGRkXz00Uce1zdj5iAHuYZrmiy3Z88ecnNz2bRpk5YYvbmXgyuvvJKBA32bIkPQNRBi4mWOHTtGbW2tU/T4RYsWYbFYWtyWESPppLOSlVzBFejQodPptDkQx+giICAAVVUJCgpCVVWqq6s5e/YsAKmpqQwePJhhw4Z57yYFAjcIMfEiqqqyZ88eSkpKsNls2iggMTGRqqoqgoODGx0hOITh5MmT2rnExESuPnk1H/ERK1lpd6m32SPRN0ZKSgqhoaEEBASg0+moqanhyJEjlJaWNml7eno6sbGxbfwNCLoyQky8iCRJWh7g+lgsFp5//nmKi4u1pFqeUFBQQK+4Xjxke4hXS15lm7qN0ZbRTdY5dOgQAOPGjXM6bzab3RXXsNlsHtkkEDSGEJN2wGAwcPXVVyPLMjabzeW1b98+vv32W5d6NpuNwsJCAGYyk3d4h2EMI5jgJvubMGECI0eOJCQkpMlyAoE3EWLSjlgsFkwmE5IkaY84qqqSnp6OLMvk5uZSXV3t9AjkeG+z2VizcQ172ctIRjbax9ixY0U4AoFfEGLSjhw7doxFixa1un4aaRzkYJNikpeXxy233NLqPgSC1iLEpB1JTU3lz3/+M1u3btXmTepvytuzZw/FxcWN1peQsNH43EZsbCw33nij9wwWCFqAEJN2Rq/XM3bsWO34v//9L4cPH/ao7iEO0YtejV4/deoUf/vb35zOTZgwgYkTJ7bKVoGgJQgx8TNz5851Oj58+DBffPEFFRUVTudLKEFB0dzqPeGOO+4gISHBK3YKBM0hxMTPrF+/nuPHjyNJEjabzclVvj6rWc1IRjpl+KvPxIkTmTBhgi9NFbSSoswvKTmwEXPFacJ6DqHn1PvI/fRP9Jp6P0FRvSg7vIOCjR+CJBE/9ibC+17gb5NbhRATPyPLMgkJCdgA7jzjAAAgAElEQVRsNtatW+e2zGEOk0ceV3BFo+2sXbuWtWvXur121113ERcX5w1zBa0gJuNKgqL7UHZ4O0kT53Ny7buo1l/9fmpOH6f/jS9hrjhN3vevCzERNE9eXh7vvfdei+ud5jR69BRSSE96ctNvbyIiIoLTp08jSRJGoxGTyST22HQSEifezpGv/q4dx1xgz85YfuxnwnoN9ZdZbUaISTvSs2dPt9HlT58+zWuvvdZovQwyAPiCL9Ch48dFPzKIQQyJHcK8efMICAhAURSf2S3wPbVlpzib9QP9rnnG36a0GiEm7YDD6cxisTjtyyktLeXNN99str6ExAVcQAYZnOQk+9nPYhZjPGXk0+c/ZTjDCSfcqc7NN9/sktRL0DGxmMo59vWLJF/2eyRd5/1Kdl7LOwFnz57l1VdfBdBiiTh2EjteBoPB4x3FEhJJJNFT6sn13a/HFG1il3UXC48v5Om0p7lkwCVaX9XV1eTk5GjH7n463tf3mI2JiSEyMrLtNy9okvrBwfO+fx1rbRX5G/9LYEQ88WNv8qNlrUfqLK7XmZmZakLlHp+1XxyUSrTpF7+0X2O2sXRHPuUmzzYANuRLviSaaMYxrvnCTaCTQE4Io29siPZfvf7/jrKARMJrf93VLDW47jjnQG3hubKAJCJqT2jn4yIC0eu8F2nfn59xZ+kjP3QoGRkZrfqld6qRSeJF83zWdml2Nokjx/ut/YdbkId83bp1WmqMSlslB3MPcmPsjfSN6gs0PRKprq6moqJCG5EUFRVpIyObCtknK8g+6ezj8iunPDeyjpY8alVW5jnFq02bdBXh4eFN1GgZ/v6MO0Mf+ZmZra7bqcSkM2Oz2Thy5Ahg94J1vBxfaneb++q/X7FiBUVFRS7tbmQj/ehHVUEV+wv2N9p/RESEJir12w8JCXHqq/5PgMrKyjbdt+OeG/Lwww/TrVs3p3Mijm3nRohJO2GxWNiyZQsHDhxoUzvx8fH079+fanM1n5/4nB2FO1gyfQmhFaH06dPHaTRiMBgIDQ0lJCTEK3FexZdd0BRCTNqJgIAAoqOj2ywmJwpO8HXB16xlLQkkMIc5nNx1kurqag4ePAjYhSQlJUV4xAraFSEmPuTUqVMeLf02x/jx45EkCYvNwvyd84k0RvK3fn9jcMRgbVWoqKiIuLg47bGpW7duWK1W9Hr3+XYEAm8jxMSHxMbGak5qzT0iqKpKVlaWFj6xqqqKVatWAbBhwwYArFippBIdOn7Y+QPHONZk1LUvvvjCY1vnzp1L3759PS4vEDREiIkPMZvNPPfcc9qxTqdz+6rvf1JSUkJiYiJBQUEEBARQW1ur1dej5y7uIpts9rCHr/iKBBJI16WTZkvTcus0x8SJExk5cqRLLh1HnNjGVoMcE7ciMbrAHUJMfIjRaNRGJllZWS4xYK1Wq+Yd60gbumPHDlJTUwHcxoutH/XeZDWxvWg7awvW8t+T/yU6OJq+5X25gAsIJLBRu5raFNgcs2fPJjIyEqvVqt1DfHy8SEEqEGLSXkiSpC0HN4bNZmPPnj1s3bq1RW0PYxhDGEJeeR4/8ROv8RrXcz1JJDVaJzw8nPDwcKdl4OYcGKurq9m8eTNgH2U57mfChAlCTARCTDoSOp2ORx991OW8Y/TiGAkArF69mrNnz2qjlMrKSnoH92aMOobNZzazuGoxt3Fbo48+ZWVllJWVeWzbhRdeSK9evcTSsKBRhJh0AuqPaoxGIwCDBg2iqqqKgwcPcvbsWaqqqoiIiCAiIoLkvGQmMpEP+ZB5zHPZBOgJjz/+OAEBAU7nsrOzvXI/gnMTISadjJdffrnR7Hy5ubl0796dhIQErtRdSVB5EIurFvN47ONEh0ZjNBopLy/XHNiCgoLo1auX0zyMYxS0YcMGbDYb8fHxDBkypD1vUdBJEWLSAVFVlZdffrlFjyEOHDmGAYYwhJOc5F8n/8Uc5qDHdb6mR48eLitKDjd/nU5HYGDjE7kCQX2EmHRAJEnioYcecjpnMpl4/vnnW9zWFKawhCV8y7dcxmUu17ds2eJy7oYbbiAtLa3FfQm6NkJMOglBQUEsWLAAVVUpKCgAnDfsHT58GJvNxtGjR1FVVUufoUPHLGbxHu+xilVMZSo6XPfpPPzww07HjuTo9Vd4qqqq3I6WVFUlLCxMeNt2cYSYdDIkSXJKX/HWW29p4tIYQQQxj3l8yqf8j/9xNVdjaPDR//Of//SKfeeffz6/+c1vvNKWoHMhxKSTsXXrVrdJzpsjhBBu4iY+5mM2sYkJtGwT4MCBAykrK2s2vsh5553XYtsE5wZCTDoZo0ePZvTo0drx9u3b+frrrz2qa8BAPPFU0Fjwo8bJysryyDYR8rHrIsSkE3P27Fl2795NQkKCFnmtIQMHDiQ1NZXTp0+zKncV+4v382zfZ4kLikOn02E0GomKinLaKxQREUGfPn201R0HIp6JoCmEmHQyDh06xOnTpwGwWq3ExsZy5swZUlJS6NHjV29Xx+RscHAwq1ev5lT1Kd7kTa7nek4dPMUpD0MwGo1G/vCHP2AwiP8qgqYR/0M6GSUlJZw6dcopxOKxY8cAu9A0RkBaACNKRvDUZU8RGBjotBO44ebDiIgIJ2ESCDxBiEknY/jw4U7H3333nUfpMsoOlLGDHazfuZ4hfe0erQ5BslgsZGRkeCW0o6DrIsSkk3LkyBH+85//eFw+hhgyyGDBngXM2TPHxdfEYrHQv39/p3OOZOqOgE1nz54lPz/fJdh1w0DUSUlJwuekCyLExI8UFRVx8uRJ9Hq9NvkJuHxZ3Z1zBDICGDp0KFarlcOHD1NVVdVofxOYwPu8z252My1mGgaDQWv3559/5ueff3bqx100/IYkJdnDHNQPsjR79myvpqgQdA6EmPiR0tJScnNzneYr6kc3cxfxrP75hIQEamtrMRgM1NTUuPXxUFWVn376CZvNhh49F3Ihu9nN8KLhLmVbw4kTJ1zO/etf/2q2XmxsLHfddZeI2nYOIcTEj6SmpmpR1dpKU8u2U6dOBSA/P5/nPniO8tryRtuJj49n/vz5TqLlSR8CgRCTToi7OQuLxUJtbW2jCb0c+2cmz5hM5tZMHrjmAa3MoUOH+OqrrwAoKChg7dq1TpHU9Hq910RPcO4ixKSD8PTTT/u8Dxs2FrOYvvTllVdeIS0tTRuB9OnTR3vUOnjwoFOMV6vVSllZGaqqEhkZ6bQ3SCBwIMSkg/DnP/8Zq9XKkSNHWLJkiU/62MQmTJgYjd0dvyUJwRyxX+uHLEhISHAa+dx0000uKT8FXQchJh0Ex+NEWloaCxYsIDMzU3v08JQePXoQFBQE4DLfUWwtZkvhFv4S9xd6GHqgqqrbydPrr7/eqW79EASZmZlasKTu3buTnp7uNCksgkp3bYSYdFAyMjLIyMjwuLy7yVGLxcKzzz4LQCX2BOTbC7YziEFamcDAQAYPHsyMGTOa7cNms4kJWEGjCDE5x8jLy2PRokWYTCan86GE8lt+y+d8Tg45PDrwUSaMmqAtS9dfona8QkJCRJY/gccIMTnHSEhI4O6779b23Lz66qvatSSSuIu7+I7veCDrAeZmzaU73bXr8fHxhIeHaw500dHRQkwEHiPE5Bxg6dKl7Nmzx6OyAQQwgxlsYhNf8iXzmKddKygocInaFhAQoKW8yM/Pp7y8XPPAHT16tNhNLNAQ/xPOAa666irS09Nd5jMqKyt56aWX3GbqG81oMsnkF34hlcZ9SDZt2kRUVBRxcXGUlZVp4hESEiI2BgqcEGJyDhMaGspTTz2lHVutVk6fPs2bb76JHj3Tmc5ylnMbtxFBhNs2evXqRe/evRkzZozwgBU0iRCTcwir1cr27duxWq18//33zZZPI41RjGIJS7iN2zBidClz4MABDhw4QG5uLlVVVWRmZrqU6devHxdeeKFX7kHQefG6mMiyPBJ4oe5wLvAWUAbMB4KA/zR2rChKy7NOCTRUVeXIkSMoiuJxnTGM4RjH2M52xjCm0XINAy9df/31WpiBiAj3oxpB18IXI5P+wHRgKrAbGI9dNOYCvYFHgcBGjl/zgT3nBIqisH79epfzDgez6upqNm3axJkzZ1rUroREP/pRSGGL6hUVFXHeeecRFhbWonqCcxevi4miKIsAZFm+AMhUFGW/LMsBwI1ArKIo+5o4FjRCr169uOwy54x8jonVrVu3uvVm9ZRKKgmlZd6ra9asYc2aNU7nxo4dyyWXXNJqOwSdG5/Mmciy3B/oCRypO6UHwgBTE8fCF7sJQkJCCAkJcXvtyiuv1MIZNHQ8q7+D+MiRI3zzzTcu9cspJ5HENtu4adMmNm3aRFJSErfffnub2xN0LiR3y4ZtQZblYGAxMA+4BdiA/TFmElABbGzsWFGUZxtrNzMzU23sy+QNTCaTtq+lK7T/2WefaSKzmMWcz/kMoOmVmtGjR5OcnNwmO5uio/2OOlr77dFHVVUVGRkZrYpY5YuRyeNADPZJ2A3Ag9hHIXcANcDrTRw3iS+XJX297NlR2l+8eDEHDx50Oted7hRQwP1T7mfYsGFO1wIDAzXfko5yD121/fbow91qnaf4Ys7kKeCpeqc+bFBkTjPHAi/SMGq91Wp1KTOCESxkIR+u/pDVq1e7befaa6/1iX2CcwfhZ3KOYLVaee6555wCTXtKNNFMZzqf8il3cRfjRo7j0ksvdSmXnZ3tDVO7HD+bfuaQ+RCVtkoSDYlMCJnA8orlTAyZSHe9fW+UyWbi26pv+U3obzptXFwhJucIer2eP/3pTwDU1tbyzjvvYLFYtElYq9XaZOT6oQwlm2x2spOAbQFs27atTfZcfPHFjB8/vk1tnCsMCxpGlD6Ko5ajjA0ey+bqzdiwOZVZX72eC4Mu7LRCAkJMzkkCAgK49957tePi4mJef/31ZusNYQg72alFYvMUnU5HamoqkZGRWqCk6upqqqurCQ4ObrH95zpjgsfwXeV32vER8xGOmo9iVs0kG5MZEjjEj9a1HiEmXYDo6GgWLFgAQHl5OeXl9uj0DdNpvL/kfb4q/Yo88uhJT+18cyt+NpvNbQjI+iEePWHmzJktKn+usL9mP5eGXkpPY0/+V/4/0gLSCJQC/W1WixFi0oGpqqriH//4R7v2eQVX8DEfcyu30oMezQqJA1mWiY6O1o49Ha7XL3fqlGfJ1M81jJJRmzsxYkSicz7qCDHpwISEhGgjiubIzs4mLi6OrVu3UlVV5TRfAu7TY6iqytGjR53akZE5y1n+p/sfz/R8RkuK3hyFhYV0797dSRwaCookSVpsFAf1xaozzxe0lPqCcX7g+ayrWoeERLwhngApoImaHRevO635iszMTNU8sOUrFZ6iP6THmuK6bNpe7Wd+lknx4WKf9d8SNrCBgxzkVm71etvBkcHIk2QkJHQGHVF9ojQR8fdn0NHbb48+jFnGDuW05jNGB7dsYrAlZOuyGRDsQ4emZtofMmsIxcXuxaS+4DveO+YytmzZ4uKE1hYyyWQnO7mFW7zWZn2qS6qpzqoG7A5xo9JHaUGW/P0ZdPT226OPTDqQ05qgdYSGhrqkisjNzeWjjz5qNxvMmPmO75jP/EaDJXlCcHAwd9xxh9v5FqPRKHYan6MIMenA9OvXr0VzJg436w0bNlBW9mtoGFVVPXKTVlGxYiWKqNYZDGKTXxdGiMk5SPfu3bVkWY7RwfTp0wF7Zr76QlOfAALQoaOaakL4dVPlkCF2v4fS0lKnQEiOCPiOnwBXX321929I0CkQYtKJefbZZ1323rSFs5zFiJFAnH0c9u7dq70PCwvDZDJhtVoZPHgwUVFR6HQ6goKC6N69u4hW34URn7yfKS4upqSkBL1e7+REtmvXLn7++ed2tSWffGKIQY++0TIVFRXa+/oiAxAZGanFS3G3JF3/p4P6xxaLhRUrVridcJ4yZUqLMhwK2h8hJn6msLCQXbt2YbVatccFh6AkJyd77OfhDVJI4Ud+5DM+YxKTiCa6+Ur1KCkpaXGfs2fPpnfv3kiSxIEDB5BlGcDFX8VodA12LehYCDHxM4MGDWLQoEHNF2yGlsS52LBhAz/88IPL+SCCmM98NrOZhSykN70Zy1jNtd4TgoKCCAsL46KLLnIRBndOaY4RmaOuLwNgCXyLEJMuSFxcHJdffjk2mw2z2eyUFiOAACYykQu5kF3s4jM+I4IIxjKWNNJcXL2Tk5MJCAjQHtHGjh3r02hsgo6LEJMugMViobKyElVV+eSTT1xSgLojkEBGM5oRjCCLLNawhj3sYRaz0PFrJr+Gj2H1N/xNnz6dUaNGee9GBB0aISZdgGPHjrFs2TLtODw8HEmSKC0tbbauHj1DGEI66SxkIQc5iIzcZJ0LL7yQwMBAampq2LhxI6qqMmTIEAICArQ9QaGhoV1qL05XQIhJFyAlJYWHHnoIsAuLyWRPCvDll19SXV3tURtGjAQT7BLUxx3uQg+4m6NpitDQUB5++GGRz7gTIcSkC6GqKrt27WL37t0trltAAac41WSS8+Z48MEHXc5JkkR4eDjQPgGZBb5DiEkXQpIkevXq1SoxWcc6xjDGbT7ihjzyyCNi/00XRIhJF2LTpk2UlZUxYsQIbb7C4RRms9ka3b9TSCHHOMZVXAXA1KlTtdUbnU6HXq+nZ8+exMbGts+NCDokQkw6IZmZmXz11Vft1t8JTtCXvgRgD9rz3XffNVn+2muvJSUlRdsfJOgaCDHphERGRhIREeHkdm42mzEajU2GWXRcq+8S7wn96Mf3fM9BDtKf/s2W/+qrr7jsssuIj48H3DurNYzKJuj8CDHphPTr189lMrMlk5cWi4UXX3yRmpoajwJGRxDB9VzPp3zKOMY1G72+qqqK//3vfx7ZAvalaneTs4LOhRCTLsh3331HTU0N4LrxrjGSSeY2bmMRi6illou4yGv2lJWV8de//lU7fvTRR4VbfSdEiEkX5LLLLuOyyy4D7DFKCgsL2blzJ4qiNFmvO925mZt5jdcYwQiCcc2JExER0aQz3MUXX+zkO1I/uHVRURHR0dFs27aN8847j8jIyFbeocAfCDHp4kRERBAREUFaWlqzZf/973/DGQgjjAoqXMRk2LBh6HQ6Tp8+jcFgQK/XYzAYtBUfvV5PeXk5Op1OExTHqpAkSVitVqxWq0ePXoKOR5NiIsvyAOC3QA0gAXrgv4qi5LaDbQIvc+TIEf7zn/+0uR0rVrcxT+rHX+nWrZsWUqF+NLaGkdliYmK4++67kSRJOK11choVE1mWZwO9gKcURbHVnTMAj8qynKkoStPrg4IOR58+fZp9DPEEG7YmAyiBfRL2z3/+c5v6EXQumhqZrFQUpRpAluXRgGMcfAJY72vDBL7hwQcf1OYo4NcJ2A8++IC8vDyP2rBhc9o53BCDwYCqqjzzzDPNRlq7++67hbPbOUKjYuIQkjoeA76od83kS6MEbcNsNvPCCy+4xIdt6PXaWnToKKWUbnRze91isWh9paamMnXqVHQ6neYxWz88Jfzq92K1+jaBlcC3eDoBGwXcA1QDKvBfn1kkcKK0tFSblKxPw7CGlZWVTo8v9913n/a+YUpQdy9HuYY/P/vsM5dwjFOYwid8wlzmEkOMW7sdbRw8eNDjJGGJiYkMHjzYo7KCjoenYpINvASIEUk7oqoqS5cudZnjcBdw2Ww2s2HDBpfr5eXlXrdrCEMwY2Yxi7mHezza/NcU3bt3R1VVSktLeeWVV9yWueCCCxg7dmyb+hH4Fk/FxAj8EfuKjgo+SEIrcMGxXOppoGZPY5N4g+EMZwtbKKKIRBLdlnGs6DTnvn/27Nkmr999991amAJBx6Wp1Zw5wC+KomxRFGVevfOTZVm+UlGUL9vFwi7Obbfd5lG55pZVVVVl8eLFVFZWAu73yzS3V6a6uprgYLtvyfETxymllO50b7S8t0ZFb775ZqvrduvWjYceekjsA2oHmpqA/VCW5XmyLM/AvoIjAYnAXkVRPm4vAwUtp6ioiOLiYjZs2EB+fr5P+iinHAMGgghqdRsDBgzAaDRqibtKSkqIjIx0mQ8C+2NOXFxc24wW+JQmH3MURXkfQJblHnXHZ9rDKIEr27dv10IfOsIu+pNwwgkjjH3sYwhDWtVGdnY2YA887Uh54dhp3HBeSIQz6Ph4NGciRKT9OHr0KB988IG/zWgWCYkruZKP+IgEElqcsKs+q1atavSawx1/zZo1LtdiYmK47bbbxCNMB8EjMZFlORUYDXZPJUVRxNKwj+jduzcLFixwe23z5s2sXr3a7bX6kc8kScJsNvvSTAASSeRiLuZjPuZmbm7U76QtWCwWYmLsy88jR45k8ODB2t6dhv4qAv/i6WrOu8ASxNKwXxkzZgxjxoxxOa+qKllZWciyrC2x1tbWupSTJIlly5ZRWFjoNZsyyKCSSt7hHWYy06PgSS2lqKgIgJUrV7Jy5co2tXXnnXdqj1IC7+KpmGwC3nHs0RF0LBwjEsdE5vr16zlx4oRTGccSrbdHLBISE5hAT3ryFV+RSCLTmEY4HXMpt6CgAL1er412BN6juV3DP2L3K4kF5siy/AuAoigXt4NtAuyBgxq6xYOrt2pZWZn2F3zcuHFauQ8//LDFYRpbQz/6cQ/3sIEN/D/+HxdxEaMY5ZJO1BsYjUaio6NdQhU0fORxXHOUM5lMbNu2jR49ejBr1izxiORlmlvNmdRehgjc89VXX1FcXKwdO+YKwDkheG1tLTt27NDKOa6FhIRoUcvc1WtYvuE5Tzf/gT1R18VczBCG8AVfYMbMeMZ7frNuuOCCCxg8eDC9e/duUzsg8vL4Gk8nYK9TFOWTuvdvKIpyj2/NEji48cYbPSrn7S/Kpk2bnBKat4QYYpjNbBay0GMxGTBgAGVlZXTrZp/ErT/i2r9/v1fEROBbmhUTWZYnADNlWc6vK5/ic6sEfmfs2LEue2EaEyxVVTGbzRQXF/POO+8A0IMemDFTRRUhNB/P1eFz0hjnn38+CQkJLbgDQXvjychkItAbmARYARHxRqDxz3/+063bvIREH/qwn/2MYITbuuPGjdPmM3Q6HcXFxURH2/1V6s8JJSYmipgnnYBmxURRlKdlWVaEC73AHQ8//DAffvghhw4dcrk2kYksZjFDGUogrh6sGzduBGD+/PkkJiaKOY1OjqdLw5myLL8I2ID/pyiK6/8cQafn6NGjbNu2DbCPCJp79GiORBJJIYXNbGYSjc/lOx6NgEYd9gQdH0/F5P+A54FA4FVghs8sEviNwsJCsrKyvNrmOMbxER8xgQlNhnp08PTTTzdbJigoiN///vfo9U3HoRW0L56KSb6iKLtkWTYCnq8VCjoVqampXHPNNU7n6s9dnDhxgsTERLfXHe/Lysq0DYkAscRqE7FhhHnFTpPJxDPPPONyfvr06YwaNcorfQhajqdiMkmW5UzsqS6ssiz/IBzXzj169OhBjx49Gr2u1+vdzmlUV1dTVlaGJEmEh4czaNAgbVL157yf0ZXoPFrRaY7k5GTmzZvXfEGBX/BUTEZgXxKuAgoURfH9LjJBh+eDDz7g6NGjTZb5gR8YxrBmH3HS09OpqalBr9dTWlqq5dWpH6c2ISGBkpISwsPDnbICCjoGnorJB3Vl3weeBeb6yiBB5+GWW25xOVdWVsaRI0dYunQpANVUE0VUs23l5OQ0W+ann34iJyeHW265RaQO7YB4KiYlQKmiKKtkWfbMJVPgEY6/uu5SYh49epTVq1c7pYdwlyoCoLKyku3bt2ttNtVfS87Xv1ZdXc2WLVs8qpuUlITVamWKaQpvlLzBCEYQQECjfdTH3cSqo4/y8nJeffVVl36nTZvG6NGjPWpf4Bs8FZPTwHRZlmOBUz60p8uxfv161q5dqx17K7dNRyGYYHrSk+1sZyyeRZd3lz/nuuuuo3//X8MbuNtLJPAvze0aHgb8HQgAyoBk4MV2sKvLMGHCBCZMmOBy/scff/Q4Kj3Y8+tERER4JEKtESqz2cyZM2cIDQ1tNChRw5WdI0eOADCe8XzMx4xmdLNpRRvjk08+cXt++PDhbvt3R0lJCQcOHGhV/9OnT2916Mg3tr/Blzlfkl+Rz+ik0VSaKymrKSM0IJQls5ZQVFnE3C/nEhEYwdMTn2ZATOd03GtuZPJX4GZFUU4ByLIcDbwJXNNkLUGbmTSpZRu228N7tCV9lJeXU1FRwfLly6EAIolkO9sZjXcfRXbu3NniOvHx8YwY4d7Fv7GRjiNWTGu4Z8Q9DIoZxLe533LLebdw8PRBLk+7nFuX3UpWURavbXuNv03+G4NiBnHv1/fy7m/ebXVf/qS535AVCJJlObneOTGN3kGpra3l1KlTbNu2jb179/rbHCeu5EoWspBwwhnIwDa19dvf/tbjkI31yxw7dozk5GQiIyPp3r3xFB2+JC0qjbSoNMxWM4dLDpMWlUZeWR7DE+wjrCBD66P9+5vmxKQUaOiSWOYjWwRtpKioiG+++Qawp9p0N+xv7AvoyRfTZDIRHBzcohgnDnrQgxu5kY/4iGCC6UvfFrfhYNGiRU7Hsixz/fXXN1vPZDLRt2/r+/UmT/74JA+OepAAfQAhxl99cM6amk5I1pFpLjiS8BDqRCQlJTF//nyftZ+dnY0syyxbtoyamhq3ZRRFabR+IolczdUsZSm/43cer+405IknnmjTY4e/eX/X+xh1Rq5IvwKAYXHD2Hx8M4NjB1Njcf977Qx03k9E4Bd0Oh1XXXWVdmyz2ZzCQl522WVYrVZt+bY+06ZN4xLrJVQcrsDY08iCic6b+rrCruFDZw/xwKoHuG7Qddyx4g7uH3U/9468l/u+uY/ymnKevOhJf5vYaqTOsgSZmZmprihf4bP2x4eMZ0PVBtG+g0PY/Z3rkE57bwn2LGd5m7f5Hb8jlFCX62qUal837AbeDCHb6S/G67UAACAASURBVD4DP/Qxs9tMMjIyWvVb71RikpGR4bP2ff1XsbO1f+TIEaqrq538OY4cOYLRaOTo0aPa+ebc6RtjJSvRo2c601tV35E/R5IkoqOjGT+++fCQne0z8EcfmZmZrRYT8ZgjcEufPn1czu3atavJOZGWMJGJvMmbpJBCGmktrh8eHk5sbCyqqhIW5p3dyIK2IcSki1NQUODkyOVupOo4ZzAYNAe7devWtanfUEK5nutZwhIu5mIyaNmoc/Pmza3q97rrriM9Pb1VdQVNI8Ski2Oz2fjxxx/90ndPejKPeSxhCYUUMo1prfaQdcecOXNISfk1/nl2drYQEh8ixKSLk5iY6HGoxIbP63l5eVRWVlJdXc2yZcta1X800dzO7XzO53zER1zDNS2OfXLvvfdqgagF/kOIiaBVfPvtt2zdutUrbQUTzI3cyBrW8A7vcD3XE0ecSzm9Xk9oaCgGg4GIiAiCg4OJi4sTQtJBEGIiaBXTpk1j2rRpbW4nJyeHTz75BB06pjCFOOL4D/9hJjMZgPOqhdVqpazM7oB95swZALKysrQNhQ1pmDr0vPPOa7O9gsYRYiLwC6qq8uqrr7rsjB7KUKKI4hM+oYACjwJRHz58mOHDhzNwoH3Pj7u0p2CPxyLwHUJMBH5BkiQeeOAB4NcAUfv376d///6cPHmSiA8j+IRPMGNmKlNd6k+dOlVbErbZbMTHxxMX5/poVJ+2pu4QNI0QE4FfsNls/Otf/3Jyxa9PGGHMZjZv8qZbMbFYLJw9e5aAgACGDRsGQFVVFaqqEhgY2Kn37nRWxG9c4Bd0Oh2PPPKI07ns7GxSUlL4+9//DsAe9pCA+/zC9dNp/Pjjj07hLGfOnCmWgP2AEBOB39m5cydbtmyhuLgYADNmfuAHssnmFm5ptN6sWbNIS0sjIKB1u48F3kWIicDnnDx50ikFaGOYMLGXvWxkI0kkcQd3NOlzEhcXJ4SkAyHEROARnqTtbCk2bBRRxGEOc4hDHOUoKaQwi1kkk9xs/TfeeKPRaxkZGcyYIbLYtidCTAQekZaW1upgzGAXjmKKya/7d5KTFFBAGGH0pjdDGcqVXOmVzH9g3/16+eWXi+j17YgQkw7IqlWr+Omnn4CWpXRoTUhGR74ei8XSCkvdY8VKEUWacOSTTwEFdKMbCXX/JjKRBBIIJthr/TrQ6XQMHTqU3Nxc7ZyqquTn52M0GunXr58QGR8gxKQD4vAurZ80vCEmk4kXX+w4WUfKKGMXu1BQKKKICCI04RjIQOKJJ4j2CZZss9nYvXs3u3fvdrmWmppKr169Wp22QtA4Qkz8SFZWlosH6OrVq/1kTeupoYZ/82+GMpSpTCWBBALpWF/WoKAgZsyYwaBBg/xtyjmLEBM/UllZyYEDB1odrayjYMBAEkkc5zhxxNGDHn4Vk5AQ53kXSZK47rrrGnWQE3gHISZ+ZMSIEYwYMUJzJ7dardhsNu3Ygar+//bOPD6q+tz/75nsIRtJgCSEsHMQkS1oABdAKVVwq2JBRESrUr3eqr8r1vZXu1h/9tWq7e2trb2tFoXeClw3RAERZN9X2Q+QBEhI2BKyJ5Nk5vz+mMxxZjJrciZked558cp8z/luh+R88l2e7/No/PWvf6WysrJN+hUREeFiQeroi9VqxWw2N5t+hWvhPKY9Rq6Wyz5tH+u19XSnOwMYQF/60oc+IZ/iJCQkuPQ3OTnZxQPbyZMnSU/3bAAnGIPhYqIoSiLwR+whRf8K/AR7rJ0ngGjgfW9pVVW7ZEwe56Dk7mzcuNElFnFb0NDQQENDQ9Dl+jd9WbFylrPkk89WtlJEEd3pTjrppJFGL3qRRpqhi6+O08QOMjIyXKxgu3XrRl1dnWHtCc0JxcjkeuBlIBE4CFyHXTTmAn2BBUCUl/RbIehPh8ZbLGJ32qsz4/z8fA4fP0xuVS55tXnkV+ezu3Y3p2tOE2YNI5lkUkklpekrlVSSSGrmcS0nJ4eePXvqZvPFxcVYrVbCwsIwm82Eh4cTGRnpsvtVUVGhp6uqqrhw4YLLNUc+k8nEmDFjPIq5EDiGi4mqqmsVRekLvAR8oqrqEUVRIoHZQE9VVQ/7SAudjP79+3uMoqdpGucqz6FeVlFLVNTLKvvO7mN73XaKK4vpm9iXId2HcG2PaxmWOoweKT3ol9gPEyYaGxvZunUrly5dMqyfFy9eJCIiArBHCMzK8m80J7gSqjWTEuAI8FhTOgyIA+p8pJsHUHEjlEfI6+rqpP6r0EYGGWTEZTA5bjJ1aXVYLBZO5p/kfP15CmsKOXHiBGsPraXAUsCVxiukkcZABnIt19KDHno9/fr1Iy4ujm7d7L9GZWVlhIeH6yMZsE/fHILhjsViwWKxoGkaBQUFVFdXB/0sHfVnYBShWDO5H9gAfAD8UlGUMdinOWeBKkVRsrFPa9zTBf7q7sgxTzp6/W3RxrFjx8jKyqKmpoY0axpxp+IYbBqMKcaEKdZEna2O/MZ89pbv5T3e4xZuIYccAI/e1jw5lJafgW/27t3b4rKhGJkcBxYBtcBY4EXso5AnAQvwZx9poYuTkpJCTEwMa9as8Xg/lVS+y3exYeMiF73WM3XqVBchEUJPKNZMjgDTnS497JbFX1roomzfvp1ly5YBdiOzpKQkzGYzZrOZxsZGfSt614VdHOEIP+SHXuuKjo6mvLwcs9lMfHx8m/S/qyN2JsJVZ/PmzS7OjsC+NnD+/HmP+fPIYxSjiMN7JL/PPvss4Pbnzp3rcZFYCA4REyEkbNu2LWRHA1JJ5ShHA8prMpmIiopC0zRsNpt+sFHTNOLj43nmmWdkS9ggREyEkDBhwgQmTJgQcP78/HwWLVoUUN544qnGvtsybdo0EhMTXdw2Ou/ggH0bOikpiTNnzjB06FAXWxQREuMQMRGuKvv27WPFihVBlQknnHrqAVi5cmWL2p09ezbp6elYLBavoTE8ndZ2x5FHREnERPBBQUGBvvagaRr19fWsW7euWb7Kykrq6+vbrF896UkVVZzmNP3o16I6/vWvfwWcNzo62qdYmEwmYmJimDRpUov60lkQMRG8kpaWxve//309nZeXx4ABA5o5Fjp37hyrV68GXJ05mUymFhl/+SOSSO7hHj7iI+Yz3+dCbDDExsbq/Xb+99hjj7kcJPRGezUmaytETLo47733Xod0gTCYwYxiFB/yIQ/zcLOzPC2hpqbGJR1oQHfBjohJF2fevHkB5zXa+tLZlcGGDRvYvHlzUOUnM5mlLGUxi5nGNHrSM+Cy11xzjW7DAt+6r3T0JyUlhYMHD/rst4MBAwaILQsiJkKI0DSNV199FZvNFrI2zJiZyUx2sIP3eZ++9OVGbqQ3vf2XbRIRb64xS0pKKCkp8esr1mQy0bNnTxETREyEEGEymXj55ZeDKuMY+Rw/fpzdu3d7XL+or68nMjKS1NRU/VDfzdU3k701m33sYylLSSaZG7mRQQzCRHMxSElJYcaMGYY8p/AtIiZCuyMjI4MJEyZw4cIFr4Zvqqq6pKOIYjzjuYEbOMxhvmr6upEbGc5wlzWVkpKSgOIA3Xbbbbqxm6ZpZGdnywjEByImQrsjISGBhIQEBg4cqBu+2Ww2Ghsb9ZAcJ0+e5NNPP21WNowwRjKSEYwgl1y2spV1rGMc4xjDmKDcR7pvg2/bto24uDiXqZGzi82GhgY+//xzGhoa0DSNBQsWdKmIgyImQofAbDYTGRnJ5cuXqampYc+ePfTo0UNf03C81I54xSZMDGr6KqKIrWxlM5vJJpsccojH/wjDEQ7DYrEAdrG4cuUKgN62o/3w8HBGjhzJsGHDCA8PJyIiossZsomYCB2KEydOcPjwYUpKSgD7S+3P41oGGTzAA5RSyg528Gf+zP3cz2AG+yznEBF3oqKiqKysdNkBGjp0KN26ddPXcboiIiZCh2LSpEm6pWllZSX5+flomobFYmHVqlU+yyaTzDSmEUccRznqV0y84S4yN998M8OGDdNHLV0VEROhw1JbW8upU6ew2WzYbDaGDh2KzWbzGRO5hhp2sYuZzDSsH5s3bw7YRuaxxx6jT58+hrXdnhAxETosK1asoLCwMKgyW9mKgkIf+pCYmEh9fT21tbVB1REbG0tOTg6RkZHk5OTo6yZt4baxPSNiIrRLLBYLv/vd7ww1equnnv3s53EeB6C8vLxF9dTU1LB+/XoAvvzyyxbVMX/+fNLS0lpUtr0iYiK0S6KiogI2ejtw4ADLly/3mWfAgAHsse1hbN1YnrjxCQCPRnHOOzSOzxkZGcTE+A8YJiMTQeiAVFdX88YbbwSc/2jeUf7BP5jJTD766KOAymRkZBAVFcVtt91G797+TfS7OiImQofj3LlzXLx4kbvvvrvZvZMnT3p0BXCRi8QSSx8CX/wsKioC4J133uH2228nJyen5Z3uAoiYCB2OkpISzp496/FeRUUFZrOZsLAwl3jJGWRgw8Z+9jOa0V7rTk5OZtiwYWiaRmZmpku8YsE3IiZChyMyMpIDBw74zOO+cBtOON/n+7zHe/SiFxlkeCxXWlrKli1bAOjTp4+ISRCImAgdjsTERHr06KFbn9psNkpLS13yOHuEKy0t5cqVK/SkJ9OZzlKW8hRPuZzTGTduHH379nVZhDWbzeTl5enp5OTkgDyuufPee1tYvfoQFy9WcsMN/cnKSiEv7xImE/z2tw9QUlLNj370PyQkRPPCC3cweHCvVvzvXD1ETIQOR3p6Ok8//bTX+6dOndIXWevq6lzuXcu15JLLWtZyJ3fq13fs2MGOHTt8tpuVlcWjjz4adH/nzbsJRUljwwaVZ5+dwmefHeDpp2/ll7/8lOPHz/P++1v4yU+moyjp/PSnH/Lmm7OCbqM9IGIidDoGDRrEj3/8YwDq6+tZu3Ytu3fv1u9PYQpv8Aa3czvhXl6BESNGkJmZ6eI0KTk5udV9i42NYtasHJYu3UVxcTmKkkZxcTkjRtgXhqOiPAdW7wiImAgdglWrVrFr165W19NIIxvYQAYZXv3GxsTEcPLkyWZm+c7C4umzI8iX+/X+/fszYIDrTtB112XyxRffcOFCBTEx3wpIebmrH9qOhIiJ0CEYMGAARUVFLr5ESkpKvJ7s9cRlLvMRHxFPPLOZ7dELGxC0eT3Ac8891yzIF3xr+LZrVz4AhYVX2L//DHfdNYqRI/tQUFDKsGG92b07n6FD06mvbwy67faCiInQIVAUBUVRWlT29OnTPPf+c3zFV0xmMmMZ61VIWsonn3yin9fxRUZGIr/73RFWrz5EbGwk11/fj6FD0/nZzz6iqsrC889PNbRfbYkpkKhl7YG9e/dq69f79lvRGrKzk9i7t0zqD1Eb9fXnsNk8x9CxWkN3dN+GjXWs4xjHeJAH6UGPFtUTEZFOZGSWzzzt/WcQCJMn9yA7O7tFStuhxCQ7Oztk9Yf6XEVHr781bWzbto1du3ZRUVEB2If+ofRa78CCheUsp4oqZjGLWGJD1ta8efO4cuUKiqIEdI6npYT657x3794Wi4lMc4SQk5WVRXR0NNXV1R4P0jnSa9euxWq1GtLmSU7yOZ8zkIHcx31ed21aS3S03VZl6dKlWK1WNmzYwDPPPEN4eNd7tbreEwttTmZmJqWlpezfv9/FAbOzQ2ZN0wwTkn3sYyMbuZu7GchAQ+p0cM0115CVlUWfPn2aHf6TU8OC0AaMGDGCESNG+MyzZMmSZiEsWsJJTjKFKYYKSU5Oji6CV65cISUlxbC6OwsiJkKLccSUsVqt2Gw2zp075xKoPDo6mu7du+t53aPnffzxx36dQQfLJS5RRBHjGd/quuLj43WXkAcOHNBHT5qmUVNTw+DBLfMh21kRMRFahGMU4YjXazabqa+vb5YvJSWl2TqJ47vZbNa9jZlMJurq6oiJidGP/oN9iuRezvHd+fOewj2sa1jHaU4ziUlk4XvnBeyLptHR0c3qtlgsREREEBMToz+bo7/OacEVEROhRcyaZfz5kZasOSxeupjXjr9GIYVMYAL3cA9RRHnMazKZCAsL04VBVVWmTu24dh3tDREToUPgyXz91KlTfH78c85xjmd4xm+0Pk3T9IiAANu3b2f79u3ce++9jBw5MjQd70KImAiGsH37dtasWdPm7Q5mMIc4xHu8xwxmkEpq0HV8+umnzUKNRkRE8OKLL3bJLd6WIv9TgiGMHz+e8eNdFz2rqqq4dOkSJpOJTZs2kZ+fH3B9AwcOxGw2Y7Vauemmm0hISPC6XrKABSw8tJBXt73K67e9zrwx8/R6nNc2uvrWbagRMRFCRlxcHHFxcQD069cP+HZX5+zZs1gsFmw2m35t48aN+u5Obm6uXk9eXl5A7c1kJs9+/iwrP1/JcIYb+zAhxPF/NHXqVK677rqr3JuWI2IiGMalS5coLy8nIiJCD97tHFjcmZiYmGZm56NHj6Z///562tNRD/dr7ul7S+5l9qrZjBo9ipyIHLZu2hpw/0eNGqW7afS0w+MI+OVpJ8dkMlFaWupif+J+etj5e1hYGKNGjXIZZTmsaTsqIiaCYZw5c4Zjx47R2NhIY2OjvsWbmmpfx3B/Cd3TFouFb775xmeeQOp5Kv4p/rXvX7xueZ0JTGA0o4nE92lesMff8edbtrU89NBDDBo0KKRtXC1ETATDGDt2LGPHjmXz5s3k5ubqPlUdOP8V9oTJZKJbt24+7zvwNUIZHD2Yl+JfIs+Sx4qyFfyx+o/kkMN4xhOB3RFRSkoKY8eO1e1GHAZpVquVhIQEoqOj9elXv379iIryvN3sTFdfkxExEQxn2LBh9OnTx+N2rqZpXL58mdLS0mbi4pgmeJoeOKcd1wIJFn4Xd5FDDmtZy3KWM4MZgD1cRqChPceOHcv06dMDytuVETERDMNisbhYrzrjvGZQUFDA8ePHdXP8UNOTnsxgBq/zOg006KOTQMnPz2fx4sXMmTNHLF99IGIiGEZZWRmbNm0KKG+fPn30FzOYLeOWEkkkqaRSTHFApvbOlJSUUFJSwiuvvOLxfmJiIgATJ06koaHB75pOsOmOgoiJYBi9evXikUceCTj/0qVLOX78eAh75IrDwC1YMfFHeXk5AJ999hmfffaZoXUbQc+ePfnhD38YcpESMREMITc3l5KSkma/sJqmUV5e7nGbt3v37owbN07Pd/r0aS5cuBCS/mlNX+WUh6R+T8ydO1f3Ceu+NezpmslkIjY2Vrc78UR7XuQVMRECYufOnaxevfpqdyNorFg5ylG2sIUwwniAB9qs7UWLFrmkn3jiCZe0pmmYTCb9u2M3y5sr1fY+/RExEVywWq288cYbzSLhuePtF9vbdWd/Jm1BAw3sZz/b2EYCCdzKrQxhiOFe6YPh73//e5u3+fDDDzNgwIA2aUvERHAhLCxMj4YHrvYbjiG2pmnU19fz29/+tln59uCg3IqVRSwikkju4z7D10j84WmqFwwLFiwgNtaz82uZ5ggdgqqqKi5evOjVzuPy5cucPXtWv+4cd3fHjh1cvnzZpUx9fT1lZaEN/eCJy1ymggqe47k2HYkMHjyY2bNnt1l77Q0RE0GnqKiI7du362l3o7Oamhr9AJ7jXkFBQdt2MgCSSaaOOmqooRveLWpbQnx8vNd7hYWF/P73v/c7EomMjOTpp58mLMxzeNKOioiJoDNkyBCGDBni9X4wQ+yGhga++OIL3Sjt0KFDhvQxECKIYChD2ctebuEWQ+uurKw0pJ5XX301qPzTpk3j+uuvN6TtUCFiIuiWqI6dBU/3wb446/BU5u9Er81mIzMzU79WWVnJ6dOnQ9B7z9zMzfyDf9CDHlxD+1xjCIaVK1eycuXKkNQdHh7Oiy++SEREcJbBzeoxqD9CB+bIkSN88skn+ovv7aSuYwvT+ZqnfN6uOcfh9ZTPZrMRFhZGfX19q83sU0nlIR5iKUspoIDbuI0wrv60ol+/fiQnJ7e4fFlZGUlJSYb15zvf+Y5hrg9ETASGDx/O8OH+nQl5muZ4OswX7GfHd1VVXaZZznmLior45z//6bePzvSmN/OZz6d8ykIWMoMZJGHci9gSTp8+3eIR2siRIxk7dqzs5ggdG6vVyieffMKyZct85nMfubTFQT5fdKMbD/Ig29nOO7zD0zxtWMzhu+66C/hW9IqLi0lPT9fve1uI9bVAq2kaiYmJKIri8f6xY8da2t2QI2IiBERYWBjf+973Av6raLFYOH/+PIWFhS6BtjwtYDperurqarp160Z1dTUXL140puOAGTM3ciMnOEEhhQzB+yJzMKxYsSKgfGPGjHFJ+zvY1x5sdVqCiIngwu7du0O20Hc1qaGG85wnk8w2b3vfvn1Bl/F3dGHmzJm6i8n2goiJ4ML111/vdQvS05qJxWLxa3rvoKioiLy8PJe1Ekd4UbPZTFlZGYmJiRw6dMjwv86HOMQQhhg2xWkNkZGR+hTQ8U/TNJ599lm/i6tiASt0WqKiogJyaQh2vx/uL8LOnTtZu3YtYN/NKSwsxGw2+6zHarUG3c8wwjDju14jycnJ0Q/vOf8LCwvj+uuvx2w2ExYW5hJutL0f5POHiIkAwPnz5ykpKXG55j46OHfunP4i+1tE9IcjT3R0NHfeeSdgH7lkZGT49RPriU2bNunm/J7oS182sxkNrVUm9o5t3VmzZtGjRw+Xe+151NAWiJgIgF0o3OPT5OXlYbVa9RfYZrOxf/9+lzyB2p00NjZ6DGzeVqSSSi211FFHDDH+C3ihtLQUgL/85S/N7iUnJ7Nt27YW1RseHs5DDz3UoSMIdtyeC4aSnZ1Ndna2nrbZbHz88cdYLBb9WlVVlYvjnkC3PmtraykuLja4x8FhwkQkkZRR1iox8UVpaakuNs3a97ODYzabdTsaX6JcXV3N7t27PdY3efJkevfu3eL+txYRE8EjZrOZGTNmsGTJElRVvdrdMYQpTOGf/JNJTGIMY9rUIlbTNO644w569eqlp73l83ZN0zQKCgro06ePx/utsaw1AhETwSezZs3SPweyJrB161Zyc3PbxEl0sIxiFGmk8SVfsp3t3MqtDGOYIQuzAwcOpL6+nptuukmP/ZOamhrw4nSg1NfXM3DgQEPrNAoRE8FQFEUhPT2dCRMm6IcCobnpvKfrFy9e9Gmu7ytdVFTEqlWr/PYvjTTmMpc88ljHOrawhdu4jUEMatXCbG5uLllZWfTs2dPQszMdCRETwVBSU1P1cKDBcuzYMRdz9GDIzMzkhhtu8Fu/s1hpmsYnxz/hZ1//jIJuBfz3nf9NRV6FV1G65ZZbMJlMjBo1yqNgHDt2rMsKCYiYCF0Yk8nEfdfcx8iokfxo8Y/I/nM23+N7DGawx/yOmEAbN24E4JlnnnE5tl9TU+PT34mvLe9At9pramr00Bq+cJzwdsT0aQtETISQkpeXx+LFi692N/xyAzfQi158wAe8yIsBraO89dZb+uf4+HgaGhrYsGGDSx5nD/T+CMS+prGxURczf2UiIiKYP39+q/2UBIqIiRA0jgDfAEePHqW6utpn/ilTpgRU78WLF+nZsydHjx7VD/o5r7uEmgtcIJ74Fq2dGOWBzWhee+01v3kef/xxQ7aURUyEoPnyyy85cOAAYHfP2BnYwx62spVHeKTVTqi7d+9OTEyMy+Kyp8/uOF/3NuKwWCwed4jmzZvn1aN9WyFiIgTNtGnTmDZtmuH1tsYc/eDBg+zYsUN/WZ3/OV5Mx4vo/KLabDYOXTnERutGHuVRkmm9rcaVK1e4cuWKx3vTp093ibMcLHl5eR7j4MTEhMYQLxhETIRWs2TJEnJzc11eYAcdwTdHJZUkkkh3unvNk5SURL9+/QDv3uLKy8v1BU9vW+H5+fnNbHB8/X+5pysrKz3a8DjnGzduXJsF3nJGxERoNdOmTeMPf/jD1e5GixnIQL7iK05wAgXPHs7Kysr0qV1rGTNmTLOX3dcZJ+d0YWEhmZmZPk3u3Q8gthUiJkKrSUhI4Be/+AWAHunvartrDIYwwpjCFDay0auYtIbhw4e7vPwNDQ2cOHGiWT5PbgjcRaasrEz3H+NtqnTy5Mmgp1Fms5nk5ORWxfIRMREMJTIykpdffjnocpqmsWzZMo+LiJqmNTutbDQDGMBSlraqjpiYGP2wpK8zNlFRUc3sP9ynRb4cbzuf73EWjda4hbDZbMTGxroc7AwWEROh3ZCYmKhbz7r/ZXXfujSZTFgsFtasWWNI22GEYaN1o6na2lqPIw53MjIyuOWWlgUHC7XPlL1797a4rIiJ0C4wmUz07ds36Bdl/PjxAec9duwYQ4cOZeHChZSVlbnYhtRTTwSejbvGjx+vm/l3796dzEzPfmTFOZIgXEUcTpk0TePSpUu6J3t/NhqB2m44U1paypkzZ8jIyKB3794uZdZeWkvWmSzw4BHSOf4y2Bc4n3rqqQ7vZtFoQiYmiqL8L/A68CugAngCiAbe95ZWVbUiVP0R2ieNjY1YLBZMJhMNDQ2sX7++zfugorKc5cxhjtc8ZrNZF61Lly7xyiuvMGnSJCZOnNiGPW3fhERMFEWZAlwPvAAsAKKAuUBfP+m3PNUndF769u1L3759Afs0wdl/CsC7775LYWFhyNrfz37WsY4HeZAMMrzm87Q7tWHDhmZncXyRkpLCv/3bv3XaEY3hYqIoSjgwD/uIY4CqqocVRYkEZgM9/aQFgYULF3L27NmQt3OCE3zN18xjHqm0zG2CMwkJCUyfPh2w+3S9GoZjV5NQjEyeAf4LuANwOKcIA+KAOh/pbv4qDmVoxLq6Oqn/KrfhqH/cuHHExcVx9OjRkLVlxcoXfMF93GeIkABUVFTwwQcf6GnHYqyzZfCgQYNc/OgGS1v8nFtKKMRkFJAFjAfGKooyBvvayFmgSlGUbOzTGvd0gb+KQ7lSHuqV+I5efzBtaJqmm3y7D+n37dvH4cOHQ9K/YDjKaEI3CAAAE1pJREFUUbrTnf70D1kbjpf+1ltv1Q3Shg4dSnx8fKvq7DJbw6qqzgNQFOXnwJPAi9hHIU8CFuDPPtJCJ6C+vp4tW7Y022FpaGjg3LlzV6lXrpzmNNfQNtu4X3/9tf7ZEXAsWCZOnMikSZMM6lFoCNlujqqqrzR9fNjtlr+00MGJiopi7ty5QZc7duwY/fr1c5kquOPvIJzzdU3TvIbYiCKKelofxycqKoo5c+YQERFBfn4+AwcO1EdjJpOJ7t27+41Q2FkQOxPBUDRN4+2339btRfzhPg1qq1PGPehBPq33oG+xWHj33XeDKvPII4/oJ5A7EyImgqGYTCaeeuqpgPK6C4e39YDGxkaOHDnitdzRo0fJzc0Nqp/d6c5+PJ/3iY6OZu7cuUREROhxgB342tY9deoUgwd/6z+2W7duHTpCX7B0nScV2oxA7Sg8RbXzNCWIjIxk9OjRXusZM2aM37Y0TdPN6TVN41TpKT7926dUN1bTzW0jsa6ujr/97W96Oisrix49euh9cxcXR9oR0c9xraXnbzoqIiZCu0LTND788MOQbgs7GMEI3uEdZjGLXvTymu/s2bMtsnvZunVrs2v//u//ftUj74UKERPhqtPY2Mjrr7/e5oHNb+VWUkjhfd7nTu5kGMMMrf/ll1/uMouvIGIitBG5ubkcPHjQ5ZpjevDNN99cjS4BMJKR9KAHi1hEf/obFtQ8IiKCt99+2+WayWQiPj6eOXPmdEqTehEToU1ISkpyMS93PvHrOJtTVFRERkZGQNu/3pwH2Ww2vvzyy6D6lkEGfemLisooRgVV1hNms5m4uDisVmuzmDkVFRUsWbLEJX8gO1g33XQTWVlZre5bKBExEdqElJQUUlJSfOaJjo72GGvYk1NmR1wdd1rqLnIsY1nOcsIJ51qubXG4i3HjxgGuJvTOHtK8Wb86i46nUYu//7v2gIiJ0K44ffo077//fpu3O5jBPMADrGQlu9jFd/kuvQk+MNWOHTuCyh8eHs6Pf/zjTrGF3PGfQOhU9OvXT3dO7Y81a9Zw4sQJSkpKDGm7L32Zz3wOcIAP+IBBDGI60716YHNn2LBhDBw4EPjW+rUzGqd5Q8RE6FB89dVXuuezUFjLmjEzhjFcy7V8yqesZjV3cVdAZY8ePeqypX3dddfR0NCApmkMGDCgU4w+fNG5n07odNx22236KVxnGhsbWbduncsaRUlJCTabDavVis1mo6KigqqqqoDaiSKKu7mbt3iLkYwki+AXPw8dOsShQ4cAeOGFF0RMBKE94c1uIzIykjvuuMNn2aNHj7pMOxxWsZqmUVVVxaZNm1zyxxDDPdzDEpYwgQnkkON1yjN16lQuXLhAWloaI0eObBfhOtsaEROhy2AymaioqMBq/dZrdHp6uj6SGTRoEGfOnGHbtm26SfyQ6iE8xmOsZS072cmN3Eg22c1ExRFy45tvvmHr1q08//zzXcpgDURMhE6IzWbj8uXL+iE9x/fKykqWLVsWdH2ppDKLWRRRxAY2sJWt3MRNjGGMx5FKVVUVv/71r5tdDzRgViD06tWL+fPntyvjNxETocOwfPlyj64NamtrKS0tDXn7GWQwm9mc4xwb2chWtnILtzCGMZjxPwrxJiDDhw/X70VHRzNt2jSvo5r2HJtHxEToMEyYMMFj+MrTp083C5HhyQDMmyVtsCOG3vTWRWUVqzjMYe7lXpJICvhZnHF2Yzlr1qx2NdoIBhETocPQo0cPj9czMzO56aab/JYPxV91q83K69te583tb/L0wKcxH/I/Qrnhhht0D2xms5n4+HgUxfiA6W2NiIkgBIDVamXFihUUFxd7NOWfwQzePPQmc5nr050B2Hekrly5gs1mw2azUVNTI2IiCF0Fs9lMeHi4RyHR0LjIRWzYiCLKb107d+7kpZdeIjIyMhRdvWqImAhCAJhMJu68807uvPNO/VqFpYIlh5fwl91/wWwys/DahVw6cMmveb+mafz+97/npZdeCnW32xQRE6FTYrVa2bNnD6tXrzasTg2NUkrJIw8VlbOcZQADeHT4o4xNHktJcQnDhw/X8x8+fNirsFgsFn71q195vHffffdx3XXXGdbvtkLERGgXfPzxx7rpeXvBho1LXOKM05cJE/3pz2hG8wAPEEUUtYdr2czmgOuNj48nIsLVPsV5F2n9+vUeA7hrmkZDQ4NuIOcJk8nE/PnziYryP90yGhET4aqxc+dOF6/zKSkpxMbGtqpOTdOoqKigoaHBJX6N2WymsbGRiIgIKioqPJZtpJEiijjDGc5ylgIKiCWWLLIYxCCmMIUkklrs68TBpEmTdIdQnvC1NZybm+sSm8cTV2stRsREuGooikJ6erqePnPmjP6SffXVVxQWFhreZm1trUu6nHIOcpCTnKSYYlJJJYssRjOau7mbeFoeytMbK1asCCjfnDlzdJcGDi5cuNBuHVKLmAhXjaSkJJKSvjX0qq6u1l0T/uAHP2D//v0UFNhDUBcVFele1DRN4/Lly61qu4QSvuZrcsnlWq7lFm6hD30C2o0xgvT0dGbPnu31fkRExFWZqrQGEROhXVFaWsqf/vSnkNXfSCNb2MJOdjKe8dzDPUTif1owcuRIwD6ayszMbDbNMJlMnDhxgiFDhvgM2uUI7NUZDwGKmAit5qOPPnIxCW+vFFLIcpbTne78kB+SSGLAZR0e9J096cfExPDCCy/owhAdHU23bt08lu8KiJgIreb+++/n/vvvDyjvO++8w7lz50Lco+bsYQ/rWc/t3M5whrd4EXXmzJlkZGTQrVs3wsLCDO5lx0bERGhTHn/8cZe0c5iKQ4cO6V7YGxoaWLRokSFtXuQiX/M1j/M4ybRu8XLp0qUAjBo1iszMTJd758+fR1GUTjmFCQQRE6FNCGQq5LyVaxQWLKxiFTdyY6uFxJkDBw64xBl2bEnbbDYRE0EIJYFMhXyd6rVYLJw7d85FcBobGwF7uIjY2Fjd/eKHH31IUXgR39R/w976vQzWBpNDjt8+jhs3zmtwcoCePXuSkZGht5mY6LrmcuzYsU7v59UXXffJhQ5FRUUFW7ZsIT8/3+N9CxZyyeWk6SS55lySbEnkxOfwauKr9I7urbtsTEpKYsqUKURERBARESHrHgYiYiK0Gk3T+M///E+vlqWhoowyVFROcIICCuhDHxRNYaJ1IknWJLBA2eUyyijTy5w9e7ZZzOO2IioqigULFnRaARMxEVqNyWTi+eef93q/oqKCNWvW+FxL0DQtqO3llazkMIcZzGCyyeb7fL/NDM4GDRrUbH3H4WM2ISEBgLFjxzazXu3siJgIhuCIV2Oz2Vy+NzY2smrVKo4fP25YW1VUsZ/9/Af/QTTRQZefPHky48ePb3bd28Kv8yE8k8nkdV2kPftnbQtETASPFBYW8u677wZVxrG74ewVXtM06uvrW90f58NrSVoSCQ0JFIcXo5i9eygbMWIE06dP19Nd/WUPNSImgkcyMzObxfxtzctYW1vL4sWLKS4ublF5d0G6i7tY1riMiUz0GnJiz5497NmzR59eaZpGWlqaoQflYmJimDZtWod1Am0kIiZCmxATE8OTTz7pM4+zWJ06dYq6ujqXUY7jsyM9s2wmP137U96seZOhDOU6rqMf/QjDdYHTcUAQoLi4uMWC5o7JZCIqKor6+npMJhPl5eWcOHGC5ORkbr75ZkPa6EiImAjtDk3TyMvLo7y8XA/h6Y17uIfJTOYwh1nHOsopZxjDGMEI+tAn5P2sq6trtjsUiKf8zoiIidDuMJlMTJ06FYD8/HyfYgKQQAITmr5KKOEIR1jGMm7ndq7lWsP7FxcXR3h4OE8//bSLx7SuviYjYiK0a/r3799s7SYQzledByAtLk2/5ullt1qt1NXVuVzzFqzLHXfXi10dEROhU+IsIr4ICwvr0m4DjKRrnkgSBMFwREwEQTAEERNBEAxBxEQQBEMQMREEwRBETARBMAQRE0EQDEHERBAEQxAxEQTBEERMBEEwBBETQRAMQcREEARDEDERBMEQREwEQTAEERNBEAxBxEQQBEMQMREEwRBETARBMAQRE0EQDEHERBAEQxAxEQTBEERMBEEwBBETQRAMQcREEARDEDERBMEQREwEQTAEERNBEAxBxEQQBEMQMREEwRBETARBMAQRE0EQDEHERBAEQxAxEQTBEERMBEEwBBETQRAMQcREEARDEDERBMEQREwEQTAEERNBEAwh3OgKFUXJAX4BnAbeA34FVABPANHA+97SqqpWGN0fQRDaBsPFBBgFLFBV9YiiKK8DC4AoYC7Q10/6rRD0RxCENiAU05zRwPOKonwKjFJV9TBwCOgPpPlJC4LQQQnFyOQFVVWrmqY7W5quhQFxQJ2PdDd/Fe/du9fgrkr97a0Nqb99tNESQiEms4G/AcnAzxVFycY+jTkLVPlIF/iqNDs72xSCvgqCYBChEJO6pilOI/A88Br2UciTgAX4s4+0IAgdFJOmaVe7D4IgdALEzkQQBEMQMREEwRBCsWbSYhRFeRhIApbiw7jNPR2osVtT/d2BnXxrWPcP7IvBrapfUZQbgN81JecC/x2C/ju38XPgJYOfIRH4IxAJ/BX4iZHP4Fb/fzU9g2H9d2rnf4HXCZHBZFP9b2Dw71BT3SE1+nSr/wPsP2NDnqHdiImiKCOAicBB7IZshhq7udWvG9Y13fNnXBeIMd1g4HZgKnAAuBn7D8NIYz3nNkLxDNcDLwOJ2P+frjP4GZzrvzsE/UdRlClN7bzgp77W1D8WGBmK/hN6o0/n+ucb+QztZpqjqupBYBFgwr9xW9DGbm71jwL+j6IonyiK0htIN6D+/1FVtQ77L9reph+QocZ6bm3cEYJnWIv9d+Il4BOjn8Gt/olG919RlHBgHva/phajf4ec6l8EfM/o/jcRaqNP5/rHGvkM7UZMmnDYkvgybgva2M2tfg27Gv8A+A32X45aI+pXFGUwkIl92BiK/ju38UAongEoAY4AI0L0DCXAYeCfIej/M9inTzYgPQT9fwb4U1P9Ifkdwm70+XhTvZNC8AyO+l8Dzhj5DO1mmuPGYaOM3TxgAh4E/g6kAFY/7QVUv6IoMdjXMx4F5oWi/25tOIwDjXyG+4EN2OfSv1QUZQz2aY4hz+BU/xJgH7DYyP5jH3FmAeOBsUb3363+HzX13cj+Q4iMPj3Un8K3f7wNeYb2KCYasBD4C6ExdtNwNaybDzQYUP9PgB7YX/bNwHMh6L9LGyF4huPYh/C12IfALxr8DM71/8Ho/quqOg9AUZSfN5UxtP9u9a8Mwf8/hN7o07n+5UY+gxitCYJgCO1tzUQQhA6KiIkgCIYgYiIIgiGImAiCYAgiJoIgGIKIiRAUTeeb2qKdBEVRotqorRhFUWa0RVudGRGTdoyiKFuaDNVQFGWToijNrBAVRfmFoii3BFnvs02HBoPtz6PAqaY2lzddKw+2/QB5BLshVUhRFGW9qqq1QB9FUa4JdXudGRGT9s1a4GZFURKAWlVVq42oVFXVP6qqusuRbhKHiQEUHaSq6vamz4MVRemPD8NHRVEWKoqSFWz/FEUZCpxSVbUx2LItwGFo9UfgoTZor9PSHi1ghW9ZDdyH3Zz5K0VRXsFu+twdmOv+simKkob9bIoF2Kqq6l8VRfkN9nMq4cA8VVUbFUX5BbBeVdVNiqLMAqZgF4e0pjb/oarq/YqiLAD2q6q6VlGUZOCcU3NXsJ+MPtDU9izsZ0nSgP8DfBfIBn6sKMpyVVXXuPcf+8ub05T+uaqqJ5vqvlNV1Tea6v0NdrcUKdhNwa/BbhnaDfgfYB128/BGIE9V1V8pivL/gAzsZ0qeBN7GfiYoHCjEbmG9ELgI9ANQVdWmKIol8B+N4I6MTNo3O7GfB7kV+0t+B1CO/SUZ4CH/j4A/qKr6MHC3oijDgKQmM/A3gBj3AqqqLsE+AvqbqqpLVVUtB8oVRekF3Nx00teBs1Pvg8DDwDdN138KlGJ/Qa9XVfVtYC/w2yYhGeql/wmqqj7oEBJFUb7T9KwOhgPbgV9jH0X8X6AM+9mRCcAc4DNVVR8BNiqKMgRIVVX1UWAVdtHSgC+wuyW4A/vZrEVNB97qvTyfECQiJu0YVVU14AJ2/xwNwDeqqv4U+1/VGi/FNC+fe2MfNXjCiuso9X3szoV2O/WlFPtfewcHsb/o55vSlU19+z3fjmAa3er11P81jpuKooQB1zQde3fwCvZT2H/APorQmtp4GfvpY2dGYRcqx3ObnT7nNo3kHL/z9c7fm9wLhCG0GBGT9s+XwFFVVVXsI4ZF2If7jpfYWTD+BPyHoiiLgS9UVT0GVDaVeRA446EM2EcCLyqK8pem9CbsJ2MXuuU7pCjKbU3lj2MfOTnq+7OiKMuA3/LtCdPlwB8VRXlFVdXjXvrvzBzgX27XHmjqew32qcpvmp7zg6b0P4F7m545qemZSxVFWYh9FLLYrT4N+6nl2YqivIH9dC7YnQAt8tAnIUDkoJ/QjKZ1Cq1pFOF+715VVT8NUbvxqqpWhqJuP+3GYJ/SrfGbWfCKiIkgCIYg0xxBEAxBxEQQBEMQMREEwRBETARBMAQRE0EQDEHERBAEQ/j/v5GPXY+UeuIAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x10990da0>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "fig_vel, ax_vel = plt.subplots()\n",
    "ax_vel.invert_yaxis()\n",
    "# plot velocity\n",
    "vel_log.plot(ax_vel, label='Velocity')\n",
    "# plot horizon\n",
    "well_cug1.plot_horizons(ax_vel)\n",
    "# plot processed velocity\n",
    "vel_log_filter_smooth.plot(ax_vel, color='g', zorder=2, label='Processed', linewidth=1)\n",
    "\n",
    "# set fig style\n",
    "ax_vel.set(ylim=(5000,0), aspect=(4600/5000)*2)\n",
    "ax_vel.legend()\n",
    "fig_vel.set_figheight(8)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Optimize for Bowers' loading equation coefficients A, B:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "a, b, err = ppp.optimize_bowers_virgin(\n",
    "    well=well_cug1, \n",
    "    vel_log=vel_log_filter_smooth, \n",
    "    obp_log='Overburden_Pressure', \n",
    "    upper='T12', \n",
    "    lower='T20',\n",
    "    pres_log='loading', \n",
    "    mode='both')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Plot optimized virgin curve:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYMAAAESCAYAAAAfXrn0AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAIABJREFUeJzt3Xd8VFX+//FXCCUI0lYWUFdXUY+hiQSlKSyC+7UiIiAgcREIvUoHgdCrCixBBOkoUn6AhV0bgqIUJaKEJBxdkSKSxQLEIMGU+f1xBxcjZAaY4U6S9/PxyCNzbu69855Jcj9zyzk3zOPxICIiBVshtwOIiIj7VAxERETFQEREVAxERAQVAxERQcVARERQMRCXGWNGG2NaBGA9CWc9/ocx5u4ArLOwMWahMWaTMWbCpa4vkIwxVxpj1hpjthhjeninnf0eJHi/P+Sd5wNjzF+908KMMeuMMbVcCS8hScVA8ovfOsxYa5dYa7cEYJ1tgSRrbWPgZmPMzQFYZ6D0A1621t4NtDHGXMlZ7wHgMcaEAcOAxsAY4Gnvz+YAofRaJAQUdjuASE7GmPLAEqAUsNZa+5wx5k5gNs4GL85au9QY8zAwGtiTY/nRQIK1dq0x5npgIXAUqA08Za39yBizCrgeOADEW2unnCNKU2CS9/GL3nVvAh6z1v5kjEmw1lb3Pt9p4H5gOrAbGGut/Ycxpj7wiLV2iHfv4i7ge6C1tTb7Et6mpsA07+OZQMQ55qkM7LfWZhpjPgIOeqePB8ZdwnNLPqRiIKFoGDDXWvu6MeZtY8w6oDzQCqdAvAQsBWKBRkBVYH4u67sT55NwY+BRY8we4DqgPrDXWtv6PMv9GfgJwFq7EcAY87tP32c9bg40sdae9M53lXf634F1xpiaQA1rbSNjzHjgIeD1nE9ojJmPsxEP865/nrX21XNkK2WtTfdm+3/eZc/+eRjwJ+C4MeY2nEK6F4ix1h727jWI/EbFQEKR4X+fej/ztosBc4FvgeLenxW21qYZY3b5WN92a22KMeY7oAiQBvwCfAw8n8tyqUBJ4KgxphewM8fPi531eM6ZQuC1zbsRrmatjTXGtAaqGGPeB0oA+8/1hNbaGB+v5YzTxphwa22Wd49jXo6fe4CfcYrGF8aYDuT+WqWA0zkDCQU5P6UmA3d6P73WwvlEGwu0Bp49a75MY0xJ4DYf6z+Vo10F2GqtrWutnZPLcluBe72PH8Q51JQBlDXG3ApUOGvetBzLrvfm/c7b/hJ4x1p7D9AbSPSR2ZdtQFNjTDjO3sf3QJgxpogxpgJOIfgPEGmMKQbUuMTnk3xOxUBCwWRjzCferwdwjtN3A7YAb1lr9wPrgM04ReFKY0whnGPfHwB9cqzP1+iL+4H2xpjNxphVxpiq55lvHtDIe7w93lq7D3gV+CcQzf829H9grd2DcyjoTW/7c5xDNpuBKfzv+P3FmgD0Aj4CFlhrf8EplB8DbwBjrLW/ep/rQ6BjjuU1QqX8TphGLZWCxnvZ6Sggy/s1zVq72dVQIi5TMRARkcCfQDbG1AH64pwk641zSV4qEINz+duS87WttamBziMiIr4F42qiJ4CeOCesdgENcTb6T+Jc1z0I5yqMc7VnByGPiIj4EIwTyG/jbNT74Jz8SwQSgBuAit4Ta+dri4iIC4KxZ1AXWIlzmOjMeC7hONdrp+fSLpHbSuPj43VyQ0TkAkVFRfnVwTAYxaAaTlf3UsAc72BYETiX0qUZY6JwDgvlbB/yteKoqKggxL14ycnJREZGuh3jd5TJP6GYCUIzlzL5JxQzxcfH+z1vMIrBYpzrwbOBZkB/nL2ALjjjt8Tl0hYRERcEvBhYa18DXjtrUs7RI6N9tEVE5DJTD2QREVExEBERFQMREUHFQEREUDEQERFUDEREBBUDERFBxUBERFAxEBERVAxERAQVAxERQcVARERQMRAREVQMREQEFQMREUHFQEREUDEQERFUDEREBBUDERFBxUBERFAxEBERVAxERAQVAxERQcVARERQMRAREVQMREQEFQMREQEKB3qFxpi2QEOcQvMXIAxIBWKACGDJ+drW2tRA5xEREd8CXgystSuAFcaYZsB04FGcjf6TwPXAIKDYedqzA51HRER8C+ZhoubAp9baRCABuAGoaK3dk0tbRERcEPA9AwBjjAH+w/828OFASSA9l3YJX+tNTk4ObNBLlJ6erkx+UCb/hWIuZfJPKGa6EEEpBsADwNvAXcaYKJzDQAeBtFzah3ytNDIyMkhxL05ycrIy+UGZ/BeKuZTJP6GYKT4+3u95g1UM6lhrnzfGfAXE4ewFdAFO+2iLiIQ0j8fDmjVr2LBhAy+99BKFCwdrM3p5BeVVWGvbeL+nAtE5fuyrLSISko4cOULPnj1Zt24drVq1Ijw83O1IAaN+BiIiPng8HhYvXkyVKlXYsmULK1asYOXKlYSFhbkdLWBUDEREcnHgwAHuv/9+nnrqKe6//36SkpJo06ZNvioEoGIgInJO2dnZxMXFUa1aNRISEnj99dd55ZVXKF++vNvRgkLFQEQkhy+//JK//e1v9OrVi7Zt25KYmMjDDz/sdqygUjEQEfHKzMxkypQp1KhRg2+//Zb33nuPefPmUaZMGbejBZ2KgYgI8MUXX1CnTh2GDRtG9+7dSUhIoEmTJm7HumxUDESkQDt9+jQjR46kdu3anDp1iq1bt/L8889TooTPQRHylfzRW0JE5CJs376djh078tVXXzFs2DBGjBhBsWLF3I7lCu0ZiEiBc/LkSfr370/9+vW54oor2LlzJ2PHji2whQC0ZyAiBcz7779PTEwMhw8fZtKkSQwYMCDfDClxKbRnICIFwvHjx4mJiaFJkyZcc8017N69myFDhqgQeOldEJF87/XXX6d79+6kpqYye/ZsunfvTqFC+ix8Nr0bIpJvHT16lDZt2vDII49Qo0YN9uzZQ8+ePVUIzkF7BiKS73g8Hl555RX69u2Lx+Nh6dKltG/fPt+NJxRIKo8ikq8cOnSIhx9+mPbt23PPPfeQlJREdHS0CoEPKgYiki9kZ2fz4osvUrVqVeLj41m7di2rVq2iQoUKbkfLE1QMRCTP++qrr7jnnnvo1q0bLVu2JCkpiUcffdTtWHmKioGI5FmZmZlMnz6dGjVqcODAAd555x0WLlxI2bJl3Y6W56gYiEietHv3burVq8fgwYPp0qULCQkJ3HvvvW7HyrNUDEQkTzl9+jT//Oc/iYqKIi0tjY8++oiZM2dSsmRJt6PlaSoGIpJnbN++nVq1ajF//nyGDBnCrl27qF+/vtux8gUVAxEJeWcPLBcREcHKlSsZP348ERERbkfLN1QMRCSkvffee1SrVo25c+cyefJkduzYQWRkpNux8h0VAxEJSceOHaNTp07ce++9/OUvf+GLL75g8ODBGlguSPSuikjIWbduHT169ODkyZPMmTOHrl27ajyhINO7KyIhIyUlhVatWtGiRQtuv/12EhMTNcLoZaI9AxFxncfjYdmyZfTr149ChQqxbNkynnjiCY0ndBkFpRgYY6oB7YHngCVAKhADROTWttamBiOPiISuAwcO0LVrV95++23atGnDzJkz+fOf/+x2rAInWHsGsUAn4BlgEFAMeBK43kd7dpDyiEiIyc7OJi4ujmHDhlG6dGlee+01mjVr5nasAivgB+KMMS2BKsAi4EZr7R4gAbgBqOijLSIFQHJyMnfffTd9+vShffv2JCUlqRC4LBh7Bu2A+4E0YLV3WjhQEkjPpV3C14qTk5MDGvRSpaenK5MflMl/oZgrkJkyMjJYsGABL7zwAhUrVmTRokXUqVOH7777ju+++86VTIESipkuRDCKwTHgMM6hnyhjTC2ccwMHgTRjTJT3Zznbh3ytONQ6miQnJyuTH5TJf6GYK1CZ4uPj6dixI3v27GHAgAHExsZyxRVXuJopkEIxU3x8vN/zBqMYzAWWAVcCDwH9cfYCugCngbhc2iKSz5w6dYrY2FimT59O1apV2bFjB7Vr13Y7luQQ8GJgrf0UaHvWpC05Zon20RaRfGLz5s3ExMRw8OBBYmNjGTJkCEWLFnU7lpyDenKISMCdOHGCbt260bhxY8qXL8+uXbsYOXKkCkEIU6czEQmoN954g27dunHixAlmzZpFjx49CA8PdzuW+KA9AxEJiKNHj9K2bVuaNWtGtWrV2LNnD71791YhyCO0ZyAil8Tj8bB8+XL69euHx+NhyZIlREdHayiJPEZ7BiJy0Q4cOMADDzzAk08+SdOmTUlOTubJJ59UIciDVAxE5IJlZ2cze/ZsqlWrxu7du1m/fj0rV66kQoUKbkeTi6RiICIX5MxQEr1796ZNmzYkJibyyCOPuB1LLpGKgYj4JSMjgwkTJlCzZk3++9//8v777zN//nzKlCnjdjQJAJ1AFhGfdu7cSadOndizZw9PP/00Y8aMueihJCQ0ac9ARM7rl19+Ydq0adSpUweAHTt2MG3aNBWCfEjFQETOadOmTdSoUYPly5czduxYdu7cqTGF8jEVAxH5nePHjxMTE8M999xDxYoVWbt2LSNGjKBIkSJuR5MgUjEQkd+sX7+eKlWq8OqrrxIXF8eHH35I5cqV3Y4ll4GKgYiQkpJCq1atePTRR6lZsyaJiYn06NGDQoW0iSgodDWRSAHm8XhYtGgRAwYMIDw8nOXLl9OuXTv1IC6Aci0GxpiiQBvgr4AH+A+w2lqbGfxoIhJM+/bto2vXrrz33ns88cQTPP/885QvX97tWOKS8+4DGmOqA9OAHdbasdbacTg3rn/eGKODiCJ5VFZWFs899xzVqlXDWsuGDRtYvny5CkEBl9sBwcJAP2utPTPBWrsH5zaWushYJA/avXs39erVY+DAgXTs2JHExEQeeOABt2NJCDjvYSJr7a4zj40xjwObgQnALmttXPCjiUignD59mvHjxzN58mRuuukmtmzZQoMGDdyOJSHE30sFHgKqADuBusGLIyKB9vHHH1OzZk0mT57M0KFD2bVrlwqB/IG/xaAk0BPYiq5AEskTfv75Z3r16sXdd99NyZIliY+PZ9y4cURERLgdTUKQvxv2QcDVwHFgcvDiiEggbNiwgW7duvHjjz8yffp0+vbtq9tPSq7OWwyMMRWARtbaVdba/+BcVgpw8LIkE5EL9v3339O3b19WrFhB06ZNefHFF7nxxhvdjiV5wHkPE1lr/wt8bIx52hijA4wiIezMfYgjIyN56623WLRoEe+8844Kgfgt13MG1trD1trngExvUbjuMuUSET+duQ9xdHQ099xzD0lJSXTo0EG9iOWC+HUC2Vq7w1sU7jDGdDfGFAtyLhHxISsri1mzZlG1atXf7kO8atUqKlas6HY0yYP8OoFsjKkL3AKE4QxLkRXMUCKSu8TERDp37sz27dvp2rUrU6ZMoXTp0m7HkjzM30tLh3q/e4BfNDaRiDtOnz5NbGwst99+Oz/++CObN29m7ty5KgRyyfy9tPRPQA/gFE5BWHq+GY0xdYDRwH5gMTAGSAVigAhgyfna1trUi3gNIgXCtm3b6Ny5M9ZaBg8ezMiRIylevLjbsSSf8LcYJAPPAul+zFsTGGStTTTGTMPpo1AMeBK43kd79gWlFykAfv75Z0aMGMHs2bOpVasWO3fupGbNmm7HknzG38NERYAhQCzOp/7c3A70N8asB2p6B7dLAG4AKvpoi8hZ/vWvf1G1alVeeuklpk6dyvbt21UIJChy63QWDfzHWrvNWvvUWdObGGOaW2vXn2fRgdbaNO/hoo+808JxhrRIz6VdwlfY5ORkX7NcVunp6crkB2Xy35lcP/30E5MmTWLDhg3UrVuX+fPnc9111/HVV1+5limUKFPg5TZq6TJjzFPGmIeAwzhXEl0NJFhrX81lne2AeUA5YJQxJgrnMNBBIC2X9iFfYSMjI/17VZdJcnKyMvlBmfyXlJREfHw8/fr1Izs7m4ULF7reZyAU3ytl8k98fLzf8+Z6zsBauwjAGFPO2/7Jj3Wmew8RZeLc+2Aizl5AF+A0EJdLW6TA2r9/P127duWjjz6iVatWzJo1S30G5LLx6wSyn0XgzLxL+f3VRtE5ZvHVFilQsrKymD17NiNGjKBEiRKsX7+eRx55xO1YUsD4dQLZGNM02EFECqKEhATq169Pv379iI6O5o033lAhEFf4ezXR7caY140xI40x1wY1kUgBkJ6ezsiRI6lVqxYnTpzgww8/5IUXXuDKK690O5oUUP6OTTQNp5PYzcAmY8zEoKYSyce2bNny253HhgwZwueff87dd9/tdiwp4Pw9TPT/gEnAAv43RpGIXIATJ07QvXt3GjZsSOnSpfnss88YP3687jwmIcHfHsgdrbUnzmoPC0YYkfzqtddeo0ePHhw/fpwZM2bQq1cv3XlMQoq/5wzuO/PAGDMnSFlE8p2UlBRatWpF8+bNqV69OomJiboFpYQkn3sGxphGwMPGmCPe+XXrJBEfPB4PCxcuZODAgYSHh7N8+XLatWunG85IyPLnMNHfcAaUa4xzH4NnghlIJK/76quv6Nq1K5s2beKJJ57g+eefp3z58m7HEsmVz2JgrR1jjLE+hqAQKfAyMjJ49tlnGTNmDBUqVOCtt97i//7v/9yOJeIXf3sgqxCI5GLnzp107tyZhIQE+vTpw7hx4yhZsqTbsUT8lmsxMMbcaa39xBjT8Ozp1toPgxtLJG84efIko0aNYsaMGVSrVo3t27dzxx13uB1L5IL52jOoB3yCc77gDA+gYiAF3jvvvEPXrl05cuQI48ePZ+DAgRQpUsTtWCIXxdeopTO938dcnjgioe+HH37g6aefZtmyZTRq1Ii3336bW265xe1YIpfEr3MGxpjCQHGgInDAWvtrUFOJhCCPx8Mrr7xCv379yMjIYP78+XTs2JFChfztriMSuvz9K14MVMYZimJh0NKIhKj9+/fzwAMP0L59exo1akRycjKdO3dWIZB8w9+/5Exr7efW2g1AdjADiYSSrKwsZsyYQdWqVdm9ezfr169nzZo1VKpUye1oIgHl79hER4wxL3gfpwQrjEgo+eKLL4iJieHTTz+le/fuTJo0idKlS7sdSyQo/O1nMMwYUwXwWGvz7h2fRfxw6tQpxo4dy/Tp07npppvYsmULd911l9uxRILK3yGsuwHTgWeNMT2DG0nEPZs2baJGjRo8++yzDB8+nM8//1yFQAoEf88ZNLDWPmCtfQCoG8xAIm746aef6NSpE/fccw/ly5dn165djBkzhmLFirkdTeSy8PecQTFjTD2cDme6E4fkGx6Ph9WrV9O7d29OnTpFXFwc3bp101VCUuD4+xffB3gYaOZ9LJLnHTp0iGbNmvH4449Tt25dkpKS6NGjhwqBFEj+nkBOAYYHOYvIZZGVlcXy5cuZNWsWJUuWZM2aNbRo0UL3GpACzddAdQ3PNV0D1UlelZCQQExMDDt27KBLly5MmTKFMmXKuB1LxHW+9gwan2OaBqqTPCc9PZ3x48czZcoUbrzxRpYuXUp0dLTbsURChq+B6jRAneR5H3zwAV26dOGbb75h2LBhDB8+nG+++cbtWCIh5UIGqrsR+AVIsdZmBjWVSAAcO3aMQYMGsWDBAurWrcvatWupWrWq27FEQpK/l5YuBcKBRcBE4ElfCxhjVgPTgDFAKhCDc1nqkvO1rbWpF5hfCiBPtocTh078blrpv5QmrJBzAvjM5aJ9+vThl19+Yfbs2XTv3l1XCYnkwt9icAw4Ya19yxjTztfMxpimwB3AQGAQUAyngFzvoz37Ql+AFDwnDp1g5l9n/m5a3/19KXN9GQ4ePEiPHj3YsGEDzZo1Iy4ujmuvvdalpCJ5R64flYwxJbwPfwTuNca8BHzvY5nCQAecT/ynrbV7gATgBqCij7bIRcnKymLmzJlUqVKF+Ph41qxZw/r161UIRPzka8/g38aYU8C7wFNAhrXW+limFzALuB84M85vOFASSM+lXQIfkpNDa4y89PR0ZfJDoDOdPHzyd+0UUmj0YCMS9ybSunVrnn76aUqVKsXevXsvW6ZACcVcyuSfUMx0IXxdTdTQGFMKaAR0Au4wxiRZa7vkslhN4Dqc+yfXNsbUwjk3cBBIM8ZE4RwWytk+5CtsZGSkHy/p8klOTlYmPwQ60/ErjgOQQQYf8AFb2cpNGRc2umgovk8QmrmUyT+hmCk+Pt7veX11OiuOMzDd3UBV4Bt89DGw1nbwLjsK6AIMxtkL6AKcBuJyaYv4VPovpan5ck36D+/P4SOHGdxjMKMnjiaiuIbNErlYvg4TvQu8B6wFhlhrPf6u2Fo71vswZ88eX22R8/rxxx8ZMGAAS5YsoUGDBmz494aQ+zQmkhf5OkykgdwlJJx9M/pff/2VuXPnEhMTo8tFRQJE/0kS8vbt28d99933u5vRd+3aVYVAJID03yQhKzMzk2nTplGtWjUSExN57bXXWLNmDVdffbXb0UTyHRUDCUk7d+7kjjvuYMiQIXTu3JmkpCSaNWvmdiyRfEvFQEJKWloa/fv3p06dOmRlZbF161ZmzZpFqVKl3I4mkq/5OxyFSNBt2LCBHj16cPToUSZMmMCAAQMoUqSI27FECgQVA3FdSkoKffv2ZdWqVTRt2pT333+fypUrux1LpEDRYSJxTXZ2NvPmzePWW29l48aNLF26lHfeeUeFQMQFKgbiiqSkJBo2bEjXrl155JFH2Lt3L9HR0boPsYhLVAzkskpPT2f06NHUrFmTlJQU3n33XZYsWcJVV13ldjSRAk3nDOSy2bx5M127dmXfvn0MHjyYZ555huLFi7sdS0TQnoFcBsePH6dTp040btyYcuXK8dlnnzFhwgQVApEQoj0DCRqPx8OKFSvo3bs3GRkZxMXF0a1bNw0jIRKCVAwk4LKzoXbtfVx1VXfeffcdmjZtyuLFi7nmmmvcjiYi56FiIAGVkZFB3brPs2tXLFdc8SfWr1/PLbfcokIgEuK0vy4Bs2PHDm67rTaffTYUiCEjI4moqEfcjiUiflAxkEuWmppKr169qFevHocPFwJ2ADPJyLiSrl3dTici/lAxkIvm8XhYu3YtkZGRLFq0iKlTp7Jr16cULXoHAEWLwosvuhxSRPyiYiAX5dChQzRv3pzHHnuMGjVqsGfPHgYOHMiNNxbm8cededq2hWuvdTeniPhHxUAuSFZWFjNnzqRKlSrs2LGDV199lX/961/ccMMNv82zeDHUqwcLF7qXU0QujIqB+G3Xrl3UqVOHfv360a5dO5KTk3n88cf/MJ5QoUKwdavzXUTyBv27ik9paWkMGDCA2rVrc+rUKT766CNefPFFypYt63Y0EQkQ9TOQXL355pv07NmTo0ePMm7cOAYOHEjRokXdjiUiAaZiIOf03Xff0bdvX9asWUPTpk3ZuHEjN910k9uxRCRIVAzkd7Kyspg7dy7Dhw8nM7MYy5Yt54kn2uk+AyL5nM4ZyG8+//xz6tevT69evahevSXp6Xs5ffoJFQKRAkDFQDh58iSDBg2idu3apKamsmHDB/zwwwKys8sxfTqkprqdUESCLeCHiYwxpYGZQFFgLjAMSAVigAhgyfna1lptdi6zDRs20LNnT1JSUhg9ejSDBw9m1qxiWOv8fO9epxfxoEHu5hSR4ArGOYM7gJFAaWA3UB1no/8kcD0wCCh2nvbsIOSRczj7BHGTJk149913ufnmmwHo0sXpMLZ3L0RGovGFRAqAgB8msta+513vUGCdtTYRSABuACpaa/fk0pYgy8rKIi4ujsjISDZv3syyZct+VwgASpeGp5+G8HAYMABKlXIxsIhcFsG6muhHIBHo6G2HAyWB9FzaJXytNDk5ObApL1F6enqeypScnExsbCwJCQk89thjDBgwgDJlyrB3794/zNugAXTseBX16v3Apb7EvPY+uSkUcymTf0Ix04UIxjmDx4DNwAog1hhTC+cw0UEgzRgThXNYKGf7kK91R0ZGBjruJUlOTs4TmdLS0oiNjWXGjBnccsstfPjhh9x9990+1zVvHkD5oGRyWyhmgtDMpUz+CcVM8fHxfs8bjD2DvcBS4BRQGxiMsxfQBTgNxOXSlgB7/fU3iI7uyenTR4mNjWXw4MHqQSwifxDwYuA9R/DgWZOic8ziqy0B8O2339KnTx/WrVsH3MuECe8zfLh6EIvIuamfQT6TmZnJjBkziIyM5KOPPqZSpZeBt1m27Cb1FxCR81IxyEd27txJmzZt6N+/P0888QQ9e+7lyJF2QNhv/QVERM5FxSAfOHHiBL179+bOO+8kMzOTrVu3MnfuXPr1K8uttzrzqL+AiORGxSCP8Hhg+HDn+/+meVi9ejWRkZEsXLiQKVOmsHr1aurVqweov4CI+E/FII946SWYOvV/t5L85ptvePDBB2ndujVRUVEkJiYyaNAgihQp8rvlOneGwYOhY8dzrFRExEtDWOcBJ07As89CVhZMm/YrBw8+y7Rp4yhXrhxr1qyhRYsW5x1ZNCwMJk68zIFFJM/RnkEeMG8e3oHjtmDt7Ywb9wxdunQhOTmZxx57TENMi8glUzHIAx577AdKl+4ENCQi4go2b/6UGTNmcOWVV7odTUTyCRWDEObxeFi0aBF33nkrp0+vISwsjlmzttOwYS23o4lIPqNiEKISExNp1KgRHTt25N577+Xrr/cydGgPOncOdzuaiORDKgYh5pdffmHYsGHUrFmTI0eO8M4777BixQquvroSEyc6J4RFRAJNxcAF5+ozAM5dx6pWrcpzzz3HiBEjSEhI4N5773UnpIgUKCoGLsjZZ+DQoUO0aNGChx56iMqVK7N7925iY2OJiIhwN6iI+DRs2DDS0tIuaJlPPvmEJUuWAPDuu++ybdu2YES7IOpncJn9vs9ABikps5g0aTQlS5bk5Zdfpm3btrpUVKQAOPN/Hip7/yoGl9n/+gxsxdpujBy5h+7duzNhwgTKlCnjdjyRfOvrr7/m+PHjF718mTJlqFy5cq7z7N27lzlz5hAeHs7QoUP5+eefmT59OoULF2bo0KF8++23LF26lEKFCnHnnXcCsG7dOkqVKkWTJk2YPXs2P/30E9nZ2dQJx/10AAAP3klEQVSqVYs777yT2NhYMjMzadKkCW3btr3o/L6oGFxmLVv+yIQJQzlx4iUiImrxr3/toHHjO9yOJZKvff/999xyyy1kZ2df9DoKFSpESkoK5cuf/+5/8+fPZ+rUqRw5coRly5bRoUMHmjdvzsaNG4mPj2fTpk3MmDGD9957jx9++OGc66hfvz61atVi2rRpREREcN999/HJJ59w//33X3R2f6gYXCbZ2dksXryYwYMHc/r0r4SFzWLWrB40bqxLRUWCrXz58nz55ZeXvGeQWyEA59BPdnY2Hu/VIa+++irlypWjQYMGZGdnExYWRqFChXK92+C111772xhjFSpUYN68eVSpUiXoRw5UDC6DhIQEunfvzscff0zbtm2ZPv1ZZs+uROfObicTKTh8HeK5FKNGjeLkyZPUq1ePoUOHEhYWxtChQ9m+fTs7d+7k+PHj1KlTh9atWzNy5EhOnjxJ7dq1fa736NGjZGVlcezYMZKSkqhSpUrQXkOYJ+f1jSEqPj7eExUV5XaM3/F1A+yzb0RfuXJl4uLiaNq0qauZ3KBM/gvFXMrkn2BkWrlyJampqWRmZlK+fHlatmx5QcvHx8cTFRXl1xUp2jMIAo/Hw9q1a+nXrx8//PADsbGxDBo0iGLFirkdTUTykMcff/yyPZf6GQTY119/zYMPPkjLli2pXr06iYmJPPPMMyoEIhLSVAwCJD09nbFjx1K1alUSEhJYs2YNGzZs4MYbb3Q7moiITzpMFABvv/02vXr1Yv/+/fTv359Ro0ZRsmRJt2OJiPhNxeASpKSkMGrUKNasWUPDhg1Zv349VatWdTuWiFygiRNhyBDnfuEXKjo6miFDhnDw4EFWrlzJsWPHyMrKYsiQIcyZM4dTp05x6tQpBg8eHNQrmi6VisFFyMjIYObMmYwaNYorr7ySpUuX0r59ew0jIZIHvfUWjBsHRYvCwIEXvnyRIkV49dVXadCgAeXKlSM6OprChQuzatUqAMaPH8+JEydYsmQJw4cPD3D6wNE5gwv0wQcfULNmTQYPHkzz5s2x1hIdHa1CIJIHZWY6hSA9HWbPhu++u/B1FC9enLS0NFJTU9m3bx+VK1fmuuuu4x//+AfgXF34888/U7x48QCnDywVAz+lpKTQvn17/va3v1GyZEk+/fRTRo4cqfGERPKwFStg61bn8YEDTkG4GC1btmTjxo2/DXeRkZHB/v37AYiNjWXx4sVER0cHIHHwqBj4kJmZycyZMzHG8O9//5t58+axbds2Qq0DnIhcuMcfhzP/ypUqQadOF76OsLAw7rrrLvbv30/FihU5fPgw+/bt4+OPPwZgzJgxTJ48mauuuiqAyQMv4OcMjDF3AlO9zSeBF4FUIAaIAJacr22tTQ10nkvx8ccf06NHDxISEujcuTMTJ04M+V+oiPivaFEYOhTatYMuXeBSzu/+/e9/58cff2Ty5Mn8+uuvjBo1itkXu6vhgmCcQL4ZuA/4O/A5cDfORv9J4HpgEFDsPO2QeOf++9//MmTIEJYsWUKtWrXYtm0bderUcTuWiARBy5awfTuMGHFxy5/Z4D/44IN/GI5i0qRJlxrvsgn4YSJr7cvW2nSgNhBvrU0EEoAbgIrW2j25tF2VmZnJP//5T4wxvPbaa8TFxfHJJ5+oEIjkc9Ong3eg0AIrKJeWGmNuBq4F9nsnhQMlgfRc2iV8rTc5OTmgOc/22WefMW7cOKy1PPbYY/Tv359y5crx5ZdfnneZ9PT0oGa6GMrkn1DMBKGZS5n8E4qZLkQwzhkUxzln8BTQwRgThXMY6CCQlkv7kK91B2OUwpSUFIYMGcLSpUu5/fbb2bZtG3Xr1vVr2YIycuKlUib/hWIuZfJPKGaKj4/3e95g7BkMA8rjFIQtQD+cvYAuwGkgLpf2ZZORkUFcXByjR48mPDycF154gZiYGMIvpguiiEgeF/BiYK0dBYw6a9KyHLPkvNj2sl98u2nTJnr37k1SUpKuEhIpoDzZHk4cOvGH6aX/UpqwQqHdifTw4cMB79FcoIaj+Pbbbxk4cCArV66kTp06fPLJJ37dbUhE8p8Th04w868z/zC97/6+lLk+9DuTBnrUgwJRDE6fPs2zzz7LhAkTKFGiBAsWLKBDhw4UKqQ+dyJyaVq1asWf/vQnihUrxvHjx7nmmmtITU1l2LBhnDp1iunTp1O4cGGGDh1KRkYGU6dOpWjRovTu3Zu33nqLpk2bcs011zBx4kR69erFpEmTKFasGN26daN48eLMmTOHtLQ0YmJiiIiIYPr06ZQuXZqyZcsG9HXk+2Lw5ptv0q9fP/bv30/Pnj0ZM2aMhpAQkYApVaoU9913H0ePHuXTTz+lbNmyFC5cmOTkZG677TaaN2/Oxo0b+eyzz6hQoQIlS5akdevWXH311edcX0pKCqtXryYsLIypU6dSrFgxrrjiChISEjhw4AAjRozg+PHjbNiwIaCvI98WA2st/fv359///jeNGzdm/fr1VKtWze1YIpLPREREAFC0aFE8Hg89evTg+++/5/jx46xcuZJy5crRoEEDsrOzuf766+nYsSOvv/46Bw8eJCwsjIyMDL799tvf1hcVFfXbIaCwsDDat29P2bJl2bNnDwcOHKBw4cIUCUKniHxXDFJTUxk3bhwzZ86kUqVKrFq1ipYtW2pUUREJirO3LeHh4QwdOpTs7GxGjRrF1VdfTXx8PMePH6dOnTocO3aMBQsWkJ6eTpMmTbj66qtZvnw511577R/WBdCmTRumTp1KdnY2vXv3pkWLFjz33HMUK1Ys8DfQ8ng8eeJr586dntxkZWV5Fi5c6KlYsaInIiLCM3r0aM/JkydzXeZSJSUlBXX9F0OZ/BOKmTye0MyVXzNlZ2V7ju0/9oev7Kxs1zIFmne76dc2Nl/sGWzdupU+ffoQHx9Pq1atmDp1Kn/961/djiUiISysUFieuGrocsnTl9McPHiQdu3a0aBBA7Kysvjggw9YtWqVCoGIyAXKk8UgLS2NZ555BmMMGzduZN68eezcuZOGDRu6HU1EJE/KU4eJMjMzWbRoEaNGjeKnn36if//+DB8+nFKlSrkdTUQkT8tTxaB69ers3buX1q1bM3nyZG64wfVRr0VE8oU8VQwqVarE0qVLueOOO9yOIiKSr+SpYrBx40b1FxARCYIwj8fjdga/xMfH542gIiIhJCoqyq9P0HmmGIiISPDkyUtLRUQksFQMRERExUBERFQMRESEPHBpqTHmz8ASIBWIsdamupwnGigDrAyFXMaYO4Gp3uaTwIshkKk0MBMoCswFhrmd6axsq4FpwJhQyGSMqQOMBvYDi0MoVzWgPfAcofF33hZoiPMB9i9AWAhkqgP0BUoAvQmN/717vVnCgE5cwO8u5IsBMMj7VQxnYzfbrSDGmBpAI2B3COW6GbgP+DvwOXA3EOFypjuAkUBpnPeqeghkwhjT1JttIKHxuwOoCQyy1iYaY6aFUK5YnI3JM6GQyVq7AlhhjGkGTAcexf2/qSeAnkANYBdOsXI7UzdrbTNjzIPAEW82vzLlhcNEFa21e4AEwNXxJ6y1u4GlOFU3JHJZa1+21qYDtYF4a21iCGR6D+dvayiwLhQyGWMKAx1wPimdDoXfndftQH9jzHqgZijkMsa0BKoAi4AbQyHTWZoDn4bC3xTwNs4Gtg/wVohkCjPGFPdmyL6QTHmhGKR7v4fj7I657UwHjpDJZYy5GbgW51ADhEAm4EcgEeeTCbifqRcwC8gGKnmnuZ0JYKC1tjMwCfibd5rbudoB9wMxwJm7rrudCWOMAf4D/OKd5HamujiHi9cA9bzT3M40DZgPjMD54Ap+ZsoLh4n2GGOicHZTD7kd5iwhkcv7KWAq8BTQIUQyPQZsBlYAscaYWji7qm7+/moC1+H809YOkUzgbHjnAeWAUaHw+wOOAYe9OaJC6L16AOfT+F0h8j5VA8YBpYA5ofA+WWu3GWNKAJ8BWRfyPuWFYrAIiANKAl1cznKGByfXHNzPNQwoj1MQtgD9QiDTXpxPJadwDl8NdjuTtbYDgDFmlDeH65m80r2HiDKB/sBE3M81F1gGXAk85M3ldiaAOtba540xXxEa24TFOB96soFmhM771An4B05h8vt90nAUIiKSJ84ZiIhIkKkYiIiIioGIiKgYiIgIKgYieYYxJtw7LMPFLHu7MaZ6oDNJ/qGriSRkGGNGAw2AFOA/1tqxxphwnE49PwKzrLW7jTE9cYYF2WStfeEin+saoK+1dvAl5B2MMwRBJtDaWvvrxa7Lz+cbCqzCGeojFWccox+AvwKfAO/hjJvV11r7dY5lw4AZOENfBDWn5E3aM5BQM95a+6S1dqy3fTVw3Frb2TscCEALa23rCykExphNZ7ettYcvpRB4PWKtfQinT8Wt53qeQPF2Lgy31u7zTroZuA2nLwdAkrW2Pc6gdz1yLm+t9eD0i2kTjHyS9+WFTmdSsAwzxnTA+ZT7KU5HnihjzEScQdOGApW97XVAMs7Imh7gqLV2pDFmEs6QE2fGIxoH3OBdZqG19j/GmOuBWGvtU94emwustW2MMQOAeJw9kf443fhftta+fo6s240xi4Fx1tqvjTHDz3qeddbaT40xm4HtQKq1dqIxpg3OsBMVgadxPpBNA9KAN621K40xC3CGOylure3ofa5InMHQzvDgDPXxOc4QKWeGSSkNpBpjbsIZOTYNWGutXWmttcaYx/3+TUiBomIgoWaStfbDMw1jzBQgwlo73DtpojGm6Zm2MaYbzlAOXwLVvePXlLHWdjDG1MTZoA4zxtQ5ax1neACstSeNMSeNMWWA+tbaZ40xrwLf4RSF+sAfioG1doAxpgrwojFmiHdj3zTH81QCFltr93rbw4E3cYrAHTgjS2YBC3HGcgIwOBvyPTme8uwbmx8EGgNfedu3GmMWAb/i9K4uhjOkci2gJc4YOr+9ZpGcdJhIQk2YH9Nyttd7N8Bv5PjZNTifwP15rtU4G9EvvW0Pzh7HSJxRH3/HGFPWGLPUWpuEM3Ll3d4fZRljzv6/+u6sQgDwszfrczjj/3wDjMUpABO950j6AceBBd69FnAKRdRZ6/kCqIAzFAJAsrX2KWttV2vtCZxB5m7DOd9S2Ju5Gs5AbyJ/oGIgoWa4MWapMWbeWdNyfpo9u/0y8LAxZhlQybvh/dkYsxRoCxzwzveRMWZ5jqtxzl7PuzjjubzsbU8C/okz2N5POUNaa48B//XuQXTwzgfOxne1MabPebLHGWNWAVNwBg8rhzO+VANgl7U2C+gGtMDZMznlfb7TOId/qnrX+QXOvSLO5wBOgekFlPAO4d0Vp+iJ/IGuJhLJQ4wxzc5z/sLXctWBX3JeZSRyhoqBiIjoMJGIiKgYiIgIKgYiIoKKgYiIoGIgIiKoGIiICPD/Aez2Vaq2t/cFAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x10990a20>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "fig_bowers, ax_bowers = plt.subplots()\n",
    "\n",
    "ppp.plot_bowers_vrigin(\n",
    "    ax=ax_bowers,\n",
    "    well=well_cug1, \n",
    "    a=a,\n",
    "    b=b,\n",
    "    vel_log=vel_log_filter_smooth, \n",
    "    obp_log='Overburden_Pressure', \n",
    "    upper='T12', \n",
    "    lower='T20',\n",
    "    pres_log='loading', \n",
    "    mode='both')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 2. determine Bowers unloading equation coefficients  $V_{max}$  and U"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "After manually select paramter U, optimze for parameter U:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "u = ppp.optimize_bowers_unloading(\n",
    "    well=well_cug1,\n",
    "    vel_log=vel_log_filter_smooth, \n",
    "    obp_log='Overburden_Pressure', \n",
    "    a=a,\n",
    "    b=b,\n",
    "    vmax=4600, \n",
    "    pres_log='unloading')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Draw unloading curve and virgin curve together with optimized parameters:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYMAAAESCAYAAAAfXrn0AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAIABJREFUeJzt3Xdc1dX/wPEX4MCJWqaZ2bA84crExJGWq+kq98xQ3DP3ptx7JOZIc5UjV5Zlmqk5UhMXy6OWW0l/KiAqiPD5/fG5+kVUuMC93gu8n48HD+65fMabi37en3HO+7gYhoEQQojMzdXRAQghhHA8SQZCCCEkGQghhJBkIIQQAkkGQgghkGQghBACSQbCwZRSI5VSn9hgO4EJXn+qlKpmg21mUUotVEptU0qNSev2bEkplUcptVYptVMp1dXyXsLPINDyva5lmR1KqRct77kopdYppco7JHjhlCQZiIzi/oAZrfVirfVOG2yzBRCita4BvKqUetUG27SV3sB3WutqQHOlVB4SfAaAoZRyAQYDNYAvgM8tP5sNONPvIpxAFkcHIERiSqmCwGIgL7BWaz1VKVURmIV5wPPXWi9RStUDRgJBidYfCQRqrdcqpV4AFgKXgQrAZ1rrXUqpVcALwBkgQGs94RGh1AbGWV7PtWx7G9BIa31NKRWotS5j2V8M8AEwGTgKfKm1/lQpVQVooLUeaLm6eAu4AjTVWsen4WOqDUyyvJ4BuD9imeLAaa31XaXULuCs5f3RwKg07FtkQJIMhDMaDMzRWm9QSv2mlFoHFASaYCaIb4AlgB/wNlAKmJ/E9ipingnXAD5WSgUBxYAqwDGtddPHrPcMcA1Aa70VQCn1wNl3gtcNgVpa65uW5Z62vP8usE4pVQ4oq7V+Wyk1GqgLbEi8Q6XUfMyDuItl+/O01iseEVterXW0JbY1lnUT/twFeAoIV0q9jplIjwG+WusLlqsGIe6TZCCckeJ/Z70HLe3swBzgPJDD8rMsWusopdShZLa3V2sdppS6CGQFooBbwG5gWhLrRQK5gctKqe7AgUQ/z57g9ex7icDiL8tBuLTW2k8p1RQoqZT6A8gFnH7UDrXWvsn8LvfEKKXctNZxliuOeYl+bgA3MJPGEaVUO5L+XUUmJ88MhDNIfJYaClS0nL2Wxzyj9QOaAlMSLHdXKZUbeD2Z7d9O1C4J7NFaV9Jaz05ivT1AHcvrjzBvNcUC+ZVSrwGFEiwblWjd9ZZ4L1rax4HNWuuaQA8gOJmYk/MXUFsp5YZ59XEFcFFKZVVKFcJMBCcBT6VUdqBsGvcnMjhJBsIZjFdK7bd8fYh5n74zsBPYpLU+DawDtmMmhTxKKVfMe987gJ6Jtpdc9cXTQGul1Hal1CqlVKnHLDcPeNtyvz1Aa/0vsAL4CmjD/w70D9FaB2HeCvrZ0j6MectmOzCB/92/T60xQHdgF7BAa30LM1HuBn4CvtBa37Hs60/AJ9H6UqFSPMBFqpaKzMbS7XQEEGf5mqS13u7QoIRwMEkGQgghbP8AWSnlDfTCfEjWA7NLXiTgi9n9bfHj2lrrSFvHI4QQInn26E3UCuiG+cDqEFAd86DfFrNfd3/MXhiPas+yQzxCCCGSYY8HyL9hHtR7Yj78CwYCgZeAwpYHa49rCyGEcAB7XBlUAlZi3ia6V8/FDbO/dnQS7VxJbTQgIEAebgghRAp5eXlZNcDQHsmgNOZQ97zAbEsxLHfMrnRRSikvzNtCidvnktuwl5eXHcJNvdDQUDw9PR0dxgMkJus4Y0zgnHFJTNZxxpgCAgKsXtYeyWARZn/weKA+0AfzKqAjZv0W/yTaQgghHMDmyUBr/SPwY4K3ElePbJNMWwghxBMmI5CFEEJIMhBCCCHJQAghBJIMhBBCIMlACCEEkgyEEEIgyUAIIQSSDIQQQiDJQAghBJIMhBBCIMlACCEEkgyEEEIgyUAIIQSSDIQQQiDJQAghBJIMhBBCIMlACCEEkgyEEEIgyUAIIQSSDIQQQiDJQAghBJIMhBBCIMlACCEEkgyEEEIgyUAIIQSSDIQQQiDJQAghBJDF1htUSrUAqmMmmucBFyAS8AXcgcWPa2utI20djxBCiOTZPBlorZcDy5VS9YHJwMeYB/22wAtAfyD7Y9qzbB2PEEKI5NnzNlFD4G+tdTAQCLwEFNZaByXRFkII4QA2vzIAUEop4CT/O8C7AbmB6CTauZLbbmhoqG0DTaPo6GiJyQoSk/WcMS6JyTrOGFNK2CUZAB8CvwFvKaW8MG8DnQWikmifS26jnp6edgo3dUJDQyUmK0hM1nPGuCQm6zhjTAEBAVYva69k4K21nqaUOgH4Y14FdARikmkLIYRTMwyD1atXs3HjRr755huyZLHXYfTJsstvobVubvkeCbRJ9OPk2kII4ZQuXbpEt27dWLduHU2aNMHNzc3RIdmMjDMQQohkGIbBokWLKFmyJDt37mT58uWsXLkSFxcXR4dmM5IMhBAiCWfOnOGDDz7gs88+44MPPiAkJITmzZtnqEQAkgyEEOKR4uPj8ff3p3Tp0gQGBrJhwwa+//57ChYs6OjQ7EKSgRBCJHL8+HHeeecdunfvTosWLQgODqZevXqODsuuJBkIIYTF3bt3mTBhAmXLluX8+fP8/vvvzJs3j3z58jk6NLuTZCCEEMCRI0fw9vZm8ODBdOnShcDAQGrVquXosJ4YSQZCiEwtJiaG4cOHU6FCBW7fvs2ePXuYNm0auXIlWxQhQ8kYoyWEECIV9u7di4+PDydOnGDw4MEMHTqU7NmzOzosh5ArAyFEpnPz5k369OlDlSpVyJkzJwcOHODLL7/MtIkA5MpACJHJ/PHHH/j6+nLhwgXGjRtH3759M0xJibSQKwMhRKYQHh6Or68vtWrV4rnnnuPo0aMMHDhQEoGFfApCiAxvw4YNdOnShcjISGbNmkWXLl1wdZVz4YTk0xBCZFiXL1+mefPmNGjQgLJlyxIUFES3bt0kETyCXBkIITIcwzD4/vvv6dWrF4ZhsGTJElq3bp3h6gnZkqRHIUSGcu7cOerVq0fr1q2pWbMmISEhtGnTRhJBMiQZCCEyhPj4eObOnUupUqUICAhg7dq1rFq1ikKFCjk6tHRBkoEQIt07ceIENWvWpHPnzjRu3JiQkBA+/vhjR4eVrkgyEEKkW3fv3mXy5MmULVuWM2fOsHnzZhYuXEj+/PkdHVq6I8lACJEuHT16lMqVKzNgwAA6duxIYGAgderUcXRY6ZYkAyFEuhITE8NXX32Fl5cXUVFR7Nq1ixkzZpA7d25Hh5auSTIQQqQbe/fupXz58syfP5+BAwdy6NAhqlSp4uiwMgRJBkIIp5ewsJy7uzsrV65k9OjRuLu7Ozq0DEOSgRDCqf3++++ULl2aOXPmMH78ePbt24enp6ejw8pwJBkIIZzS9evXad++PXXq1OH555/nyJEjDBgwQArL2Yl8qkIIp7Nu3Tq6du3KzZs3mT17Np06dZJ6QnYmn64QwmmEhYXRpEkTPvnkE9544w2Cg4OlwugTIlcGQgiHMwyDpUuX0rt3b1xdXVm6dCmtWrWSekJPkF2SgVKqNNAamAosBiIBX8A9qbbWOtIe8QghnNeZM2fo1KkTv/32G82bN2fGjBk888wzjg4r07HXlYEf0B4YBvQHsgNtgReSac+yUzxCCCcTHx+Pv78/gwcPxsPDgx9//JH69es7OqxMy+Y34pRSjYGSwLfAy1rrICAQeAkonExbCJEJhIaGUq1aNXr27Enr1q0JCQmRROBg9rgyaAl8AEQBP1jecwNyA9FJtHMlt+HQ0FCbBppW0dHREpMVJCbrOWNctowpNjaWBQsW8PXXX1O4cGG+/fZbvL29uXjxIhcvXnRITLbijDGlhD2SwXXgAuatHy+lVHnMZwNngSillJflZ4nb55LbsLMNNAkNDZWYrCAxWc8Z47JVTAEBAfj4+BAUFETfvn3x8/MjZ86cDo3JlpwxpoCAAKuXtUcymAMsBfIAdYE+mFcBHYEYwD+JthAig7l9+zZ+fn5MnjyZUqVKsW/fPipUqODosEQiNk8GWuu/gRYJ3tqZaJE2ybSFEBnE9u3b8fX15ezZs/j5+TFw4ECyZcvm6LDEI8hIDiGEzUVERNC5c2dq1KhBwYIFOXToEMOHD5dE4MRk0JkQwqZ++uknOnfuTEREBDNnzqRr1664ubk5OiyRDLkyEELYxOXLl2nRogX169endOnSBAUF0aNHD0kE6YRcGQgh0sQwDJYtW0bv3r0xDIPFixfTpk0bKSWRzsiVgRAi1c6cOcOHH35I27ZtqV27NqGhobRt21YSQTokyUAIkWLx8fHMmjWL0qVLc/ToUdavX8/KlSspVKiQo0MTqSTJQAiRIvdKSfTo0YPmzZsTHBxMgwYNHB2WSCNJBkIIq8TGxjJmzBjKlSvHf//9xx9//MH8+fPJly+fo0MTNiAPkIUQyTpw4ADt27cnKCiIzz//nC+++CLVpSSEc5IrAyHEY926dYtJkybh7e0NwL59+5g0aZIkggxIkoEQ4pG2bdtG2bJlWbZsGV9++SUHDhyQmkIZmCQDIcQDwsPD8fX1pWbNmhQuXJi1a9cydOhQsmbN6ujQhB1JMhBC3Ld+/XpKlizJihUr8Pf3588//6R48eKODks8AZIMhBCEhYXRpEkTPv74Y8qVK0dwcDBdu3bF1VUOEZmF9CYSIhMzDINvv/2Wvn374ubmxrJly2jZsqWMIM6EkkwGSqlsQHPgRcAATgI/aK3v2j80IYQ9/fvvv3Tq1Inff/+dVq1aMW3aNAoWLOjosISDPPYaUClVBpgE7NNaf6m1HoU5cf00pZTcRBQinYqLi2Pq1KmULl0arTUbN25k2bJlkggyuaRuCGYBemut9b03tNZBmNNYSidjIdKho0ePUrlyZfr164ePjw/BwcF8+OGHjg5LOIHH3ibSWh+691op1QzYDowBDmmt/e0fmhDCVmJiYhg9ejTjx4/nlVdeYefOnVStWtXRYQknYm1XgbpASeAAUMl+4QghbG337t2UK1eO8ePHM2jQIA4dOiSJQDzE2mSQG+gG7EF6IAmRLty4cYPu3btTrVo1cufOTUBAAKNGjcLd3d3RoQknZO2BvT9QBAgHxtsvHCGELWzcuJHOnTtz9epVJk+eTK9evWT6SZGkxyYDpVQh4G2t9Sqt9UnMbqUAZ59IZEKIFLty5Qq9evVi+fLl1K5dm7lz5/Lyyy87OiyRDjz2NpHW+j9gt1Lqc6WU3GAUwondm4fY09OTTZs28e2337J582ZJBMJqST4z0Fpf0FpPBe5akkKxJxSXEMJK9+YhbtOmDTVr1iQkJIR27drJKGKRIlY9QNZa77MkhTeVUl2UUtntHJcQIhlxcXHMnDmTUqVK3Z+HeNWqVRQuXNjRoYl0yKoHyEqpSkAJwAWzLEWcPYMSQiQtODiYDh06sHfvXjp16sSECRPw8PBwdFgiHbO2a+kgy3cDuCW1iYRwjJiYGPz8/HjjjTe4evUq27dvZ86cOZIIRJpZ27X0KaArcBszISx53IJKKW9gJHAaWAR8AUQCvoA7sPhxba11ZCp+ByEyhb/++osOHTqgtWbAgAEMHz6cHDlyODoskUFYmwxCgSlAtBXLlgP6a62DlVKTMMcoZAfaAi8k056VouiFyARu3LjB0KFDmTVrFuXLl+fAgQOUK1fO0WGJDMba20RZgYGAH+ZZf1LeAPoopdYD5SzF7QKBl4DCybSFEAn88ssvlCpVim+++YaJEyeyd+9eSQTCLpIadNYGOKm1/ktr/VmC92sppRpqrdc/ZtV+Wusoy+2iXZb33DBLWkQn0c6VXLChoaHJLfJERUdHS0xWkJisdy+ua9euMW7cODZu3EilSpWYP38+xYoV48SJEw6LyZlITLaXVNXSpUqpz5RSdYELmD2JigCBWusVSWyzJTAPKACMUEp5Yd4GOgtEJdE+l1ywnp6e1v1WT0hoaKjEZAWJyXohISEEBATQu3dv4uPjWbhwocPHDDjjZyUxWScgIMDqZZN8ZqC1/hZAKVXA0r5mxTajLbeI7mLOfTAW8yqgIxAD+CfRFiLTOn36NJ06dWLXrl00adKEmTNnypgB8cRY9QDZyiRwb9klPNjbqE2iRZJrC5GpxMXFMWvWLIYOHUquXLlYv349DRo0cHRYIpOx6gGyUqq2vQMRIjMKDAykSpUq9O7dmzZt2vDTTz9JIhAOYW1vojeUUhuUUsOVUkXtGpEQmUB0dDTDhw+nfPnyRERE8Oeff/L111+TJ08eR4cmMilraxNNwhwk9iqwTSk11q5RCZGB7dy58/7MYwMHDuTw4cNUq1bN0WGJTM7a20RrgHHAAv5Xo0gIkQIRERF06dKF6tWr4+HhwcGDBxk9erTMPCacgrUjkH201hEJ2oPtEYwQGdWPP/5I165dCQ8PZ/r06XTv3l1mHhNOxdpnBu/fe6GUmm2nWITIcMLCwmjSpAkNGzakTJkyBAcHyxSUwikle2WglHobqKeUumRZXqZOEiIZhmGwcOFC+vXrh5ubG8uWLaNly5Yy4YxwWtbcJnoHs6BcDcx5DIbZMyAh0rsTJ07QqVMntm3bRqtWrZg2bRoFCxZ0dFhCJCnZZKC1/kIppZMpQSFEphcREUufPmv5/vuOFC6cn02bNvHee+85OiwhrGLtCGRJBEI8xunTMGJEGN9/705cXBNatz7D1193JXfu3I4OTQirJZkMlFIVtdb7lVLVE76vtf7TvmEJ4dwMA7Zvh6lT77JxoyuGkYOnn/6JOXPK0qjRAEeHJ0SKJXdlUBnYj/m84B4DkGQgMqVbt+C772DmTAgKgqxZT+Pm9hVDhxZj6NCeZM2a1dEhCpEqyVUtnWH5/sWTCUcI53T2LPj7w/z5EB5uUKTIYWAglSvHMn/+XEqUKOHoEIVIE6ueGSilsgA5gMLAGa31HbtGJYQTMAzYvRtmzIB16yBnToMqVTT797clKuo48+dPxsfHB1dXa4frCOG8rP1XvAgojlmKYqHdohHCCcTEwJIlUKECVKsGR47AsGFX8fZuxG+/eVKzZjFCQ0Pp0KGDJAKRYVhbjuKu1vowcFgp1cyeAQnhKJcvw5w5MHs2/Pcf1KkDP/4Yx8mTXzF8+FDy5csncw2IDMvaZHBJKfW15XWYvYIRwhECA2H6dPPBsIsLtGkDPXtCXNwRfH19+fvvv+nSpQvjxo3Dw8PD0eEKYRfWjjMYrJQqCRha6/Q747MQFvHxsHEjTJsGW7fCs7kiGDHCg44dIVeu23z55ZdMnjyZV155hZ07d/LWW285OmQh7MraEtadgcnAFKVUN/uGJIT93Lpl3gqqV+9l6taF69dhWZ8ATscWZUjpDQQGbqNs2bJMmTKFIUOGcPjwYUkEIlOw9ulXVa31h1rrD4FK9gxICHu4dAmGDYPnn4euXeGll+6wYwcc+NugVfAQou5E0b5rV2rWrEnBggU5dOgQX3zxBdmzZ3d06EI8EdY+M8iulKqMOeBMZuIQ6cbRozB1Knz/PWTLBj4+5vOA2NjzeHp6Yqz/kVVbttADuH3hAv6dOtF59mzpJSQyHWuTQU/LFwm+C+GUDAO2bIHJk83vRYrAqFHQsSPkz28uExoK586do2vPnvxsGNT38MD/+ecpeucOSCIQmZC1D5DDgCF2jkWINLlzB5YvhylTzB5Cr79ujhdo1sy8KrgnLi6OZcuWMXPmTHLnzs3q1av55JNPZK4BkaklV6iu+qPel0J1wpmEh8PcuWa9oIsX4b33zFtDtWqZXUUTCgwMxNfXl3379tGxY0cmTJhAvnz5HBO4EE4kuevhGo/4esfOMQnxaIYB/frdb547B337mg+Fhw+H2rXNZwSbNpmvEyaC6Ohohg0bRvny5bl+/TpLlixh7ty5kgiEsEiuUJ0UqBPOY8MGmDWLwGIfMelADZYvh5w5zd5BPXvCc889erUdO3bQsWNHTp06xeDBgxkyZAinTp16srEL4eRSUqjuZeAWEKa1vmvXqIRIxIg32DF6FxNj1vJrrxoULQrjx4OvL+TN++h1rl+/Tv/+/VmwYAGVKlVi7dq1lCpV6skGLkQ6YW1voiWAG/AtMBZom9wKSqkfgEnAF0Ak4IvZLXXx49pa68gUxi8yuLg4WL8eJgwO5+8TkyhFEItc2lFv6Ae4fvAe8dch/Dp4PO+Bi6t5X8gwDH744Qd69uzJrVu3mDVrFl26dJHuokIkwdpkcB2I0FpvUkq1TG5hpVRt4E2gH9AfyI6ZQF5Ipj0rpb+AyJhiog2W1v+BiaebcuIEVC90mZ+L+/Fh3l3ExcYR2OUIGzh2f/lep3uR74V8nD17lq5du7Jx40bq16+Pv78/RYsWdeBvIkT6kOSpklIql+XlVaCOUuob4Eoy62QB2mGe8cdorYOAQOAloHAybZHJRUbCpEnwUpEYfLc0pWSBS/z1F+wIU3x0cgYuBwOI+nk7G2j4wHpxcXHMmDGDkiVLEhAQwOrVq1m/fr0kAiGslNyVwa9KqdvAFuAzIFZrrZNZpzswE/gAeNbynhuQG4hOop2LZISGOleNvOjoaInJCtbEdPWqG0uXFmD58vzcvu1KswJbGUo/ns/yFOc85pNw9ZsXbj6wbhhhvP3R2wQfC6Zp06Z8/vnn5M2bl2PHjvE4zvg5gXPGJTFZxxljSonkehNVV0rlBd4G2gNvKqVCtNYdk1itHFAMc/7kCkqp8pjPBs4CUUopL8zbQonb55IL1tPT04pf6ckJDQ2VmKyQVExnzpgjhRcsMNu+vtC3zGaKdW4AxMFeNzxPnID69e+vE54zHIBYYtnBDvawh1diU1Zd1Bk/J3DOuCQm6zhjTAEBAVYvm9ygsxyYhemqAaWAU0CSA8601u0s644AOgIDMK8COgIxgH8SbZFJaG32Blq2DHLnNocP9OwJTz8N+Kwwhw/fs379A8nA43kPyn1Xjj5D+nDh0gUGdB3AyLEjcc8hZbOESK3kbhNtAX4H1gIDtdaGtRvWWn9pedkm0Y+Sa4sM7PBhGDsWVq+GZ54xX3fuDHnyJFho4eNnVr169Sp9+/Zl8eLFVK1alY2/bnS6szEh0qPkbhNJIXdhE3/9BWPGmBPKvPAC+PvDZ5+Bu5Un84Zh8P3339O7d2/u3LnDnDlz8PX1le6iQtiItV1LhUgxw4Dt22Hw4GLs2wclSsCiRdCyJWTNav12/v33X7p06cLmzZtp1KgRM2fOpEiRIvYKW4hMSU6rhM0ZBvz2G1SrBjVrwvXrbqxcCSEh8Omn1ieCu3fvMmnSJEqXLk1wcDA//vgjq1evlkQghB1IMhA2Yxjw88/g7Q3vvw/R0eaz37VrT9G0Kbi5Wb+tAwcO8OabbzJw4EA6dOhASEgI9RM8RBZC2JYkA5Fm8fHmQd/LC+rVgyxZ4Jdf4O+/oUGDlM0VExUVRZ8+ffD29iYuLo49e/Ywc+ZM8j6uAJEQwiYkGYhUi4+HtWuhfHn4+GOzi+jvv8Pu3fDBBw/PJZCcjRs3UqpUKebMmcOYMWMICAigUiWZcluIJ0GSgUix+HhYswbKlYNGjSBfPti2Df7889ETyiQnLCyMZs2aUbduXUqUKEFQUBCDBg0ia0qeMgsh0kSSgbDavSuBcuWgcWN46ikzCWzfDu+8k5rtxTNv3jxee+01tm7dypIlS9i8eTPFixe3dehCiGRIMhDJMgzzmUD58uaVQIECZhLYti11SQAgJCSE6tWr06lTJxo0aMCxY8do06aNzEMshINIMhCPZRjmILEKFcxnAnnzwh9/pP5KAMxiXiNHjqRcuXKEhYWxZcsWFi9ezNNPP23L0IUQKSTJQDzEMGDLFqhcGerWhezZzfaOHVCjRuq3u337dl5//XXGjh1L//79CQwMpHbt2rYLXAiRapIMxAN27oS334Z33zVnGfv1V7N3UOIJ5lMiPDyc9u3bU6NGDQoUKMDBgwcZM2YMOXLksG3wQohUk3IUAjDHBAwfbo4cLlPmf4VC03IL3zAMli9fTo8ePYiNjcXf35/OnTtLPSEhnJAkg0wuONhMAuvWmbWDVqyAJk1SNlAssfh4qFDhX55+ugtbtmymdu3aLFq0iOeee852gQshbEqSQSb177/g52fOJ/D882bV6DZtzNHDaREbG0ulStM4dMiPnDmfYv369ZQoUUISgRBOTq7XM5mwMOjeHV57zbwlNGMGHD9ulpNOayLYt28fr79egYMHBwG+xMaG4OXVwCZxCyHsS5JBJhERAcOGQfHisHQpjBwJ//wDPXqYvYXSIjIyku7du1O5cmUuXHAF9gEziI3NQ6dOtoheCGFvkgwyuJgYmDYNXn4ZpkyBbt3MW0RDh5q1hNLCMAzWrl2Lp6cn3377LRMnTuTQob/Jlu1NALJlg7lzbfBLCCHsTpJBBhUfb14BKGXOL/zxx3DiBEycaJaRSKtz587RsGFDGjVqRNmyZQkKCqJfv368/HIWmjUzl2nRAooWTfu+hBD2Jw+QM6AtW2DAAHO+4fr1zXLSJUvaZttxcXHMmjWLYcOGkStXLlasWEHTpk0fKCOxaBGcPJnkVMZCCCcjySADOXIEunV7nt27oVIlcwDZWzacxfrQoUP4+voSEBBAx44dGT9+PPnz539oOVdX2LPHdvsVQtif3CbKAM6fN3sDvfEGnDuXjdWrzYOxrRJBVFQUffv2pUKFCty+fZtdu3Yxd+7cRyYCIUT6JFcG6diNG+YzgClTIGdOs5to9er/8Prrnjbbx88//0y3bt24fPkyo0aNol+/fmTLls1m2xdCOAe5MkiH4uLgm2/g1Vdh0iSze+jJk+Z3Wx2nL168SJMmTahXrx4lSpQgMDCQIUOGSCIQIoOSZJDObN1qzivg62uWkT52DCZMMGcbs4W4uDj8/f3x9PTkl192sHTpMjZv3swrr7ximx0Ikc4NHjyYqKioFK2zf/9+Fi9eDMCWLVv466+/7BFamkgySCdOnDAnl69dG3LkMJ8JrFgBL75ou30cPnyYKlWq0L17d8qUaUx09DFiYlrJhDNC2MC9/0d16tShcuXKDo7mYfLMwMlFRMDo0ebzgGeeMWs3+w6LAAAaGklEQVQJtWiRtkJyid28eRM/Pz+mTZvGq6++ysaNO/j88+rEx8PkyWbhurx5bbc/Ieztn3/+ITw8PNXr58uXz6rpV48dO8bs2bNxc3OjUaNGnDx5ksmTJ5MlSxYGDRrE+fPnWbJkCa6urlSsWBGAdevWkTdvXmrVqsWsWbO4du0a8fHxlC9fnooVK+Ln58fdu3epVasWLVq0SPXvkFI2TwZKKQ9gBpANmAMMBiIBX8AdWPy4ttY60tbxpFfx8WZ//cGDITLS/D5gAOTKZdv9bNy4kW7duhEWFsbIkSMZMGAAM2dmR2vz58eOmaOI+/e37X6FsJcrV65QokQJ4uPjU70NV1dXwsLCKFiwYJLLzZ8/n4kTJ3Lp0iXmzZtH3759adiwIVu3biUgIIBt27Yxffp0fv/9d/7v//7vkduoUqUK5cuXZ9KkSbi7u/P++++zf/9+Pvjgg1THnxr2uDJ4ExgOeABHgTKYB/22wAtAfyD7Y9qz7BBPurN3L/Tsac4x0KyZ2WOoWDHb7uPixYv06tWL1atXU6tWLbZs2cKrr74KQMeO5oCxY8fA0xOpLyTSlYIFC3L8+PE0XxkklwjAvPUTHx+PYRgArFixggIFClC1alXi4+NxcXHB1dU1yY4XRYsWJWvWrAAUKlSIefPmUbJkSfLZ6kGglWyeDLTWvyulXgAGAeu01sFKqWxAS+AZrXVQEu1M7b//YOBAWLwYypY15xp++23b7iMuLo45c+bc7xm0dOlSWrV68LmAhwd8/jl06QJ9+8otIpH+WHOLJy1GjBhBlixZKF26NIMGDcLFxYVGjRpx9epVDhw4QHh4ON7e3jRt2pThw4dz8+ZNKlSokOx2L1++TFxcHNevXyckJISStiodYAWXexnNlpRSuYEegI/W+lWlVA5gOuCqtfZ9THua1rrz47YZEBBg5MyZ0+axpkV0dDTu7u5p3s7du7B8eX6++qogrq7Qs+cVmja9nqqS0knFFBoaip+fH4GBgTRq1Ii+ffs+9uzDMGD69Kfp3fv/0jTbWXIxOYozxgTOGZfEZB1bxPTbb79x8+ZN4uLiyJ8/f5rnCL916xZeXl5W/Q+2eTJQSjUCtgN5AA1UxrxNVAOIAnZh3hZ6qK21HvO47QYEBBheXl42jTWtQkND8fRM2wCvnTvNSqJBQdChA4wZA1ZcnaYopqioKPz8/Jg+fTolSpRg7ty5VKtWLU1xpzUmR3PGmMA545KYrOOMMQUEBFidDOzxzOAYsAS4DVQABgC5gY5ADOCfRDvTuHzZfCC8eDFUqGA+J7B0NrCpDRt+ok2bbsTEXMbPz48BAwbIwDEhxEPs8cwgGPgowVttEi2SXDtDi4+HefPM3kEAX39tDiBzc7Ptfs6fP0/Pnj1Zt24dUIcxY/5gyBAZOCaEeDQZdPYEHTkCVaqYD2YbNACtoXNn2yaCu3fvMn36dDw9Pdm1azfPPvsd8BtLl75CpHTcFUI8hiSDJ+DmTbOfvpeXOYhs+3ZzDMEzz9h2PwcOHKB58+b06dOHVq1a0a3bMS5dagm43B8vIISwr3Xr1rF169YUreMM5SpkBLKdbdpknv2HhYGfn5kU0jrncGIREREMGzYMf39/Xn31Vfbs2UPlypWJiDBLVsh4ASGcX8JyFY4gycBOrlyB3r3h+++hRg1z9jHLmK5UMQxz3uIxY7jf1dMwDFavXk2vXr2IiIhgwoQJvPfee5QtWxaQ8QIi8/rnH0jDmDPy5YPHDVUYPHgwQ4cO5fz58yxevJjChQtz7tw5vL29H+pNlLBcxaBBg7hx44bTlquQZGBjhgHffWcmgrg4WLDAnHgmrX31v/nGHIlcvDi0bw+nTp2iW7du/Prrr9StW5evvvqKF198kdDQ0AfW69ABTp0CH5+07V+I9OLKFShRwuyskVqurubVvDXdvN977z1iY2NZs2YNTZo0eeBnCctVLF26lHbt2jltuQpJBjZ0/rx5K+aXX8zibjNnQuHCad9uRIQ5gU1cHEyadIezZ6cwadIoChQowOrVq/nkk08eW1nUxQXGjk17DEKkFwULwvHjab8yeFwicHFx4c6dO5w/fx4wy0nce/2oZdNLuQpJBjZgGOaZe79+ZnnptWvh449tt/1587AUjtuJ1p0ZNeoYPXv2YNSoUeTJk8d2OxIig7BnNQpvb2/8/f2tKu3u4+Nzv1zFoEGD2Lt3r/OWqzAMI118HThwwHA2ISEhxtmzhlGnjmGAYXz6qWFcvWr7/fzzzxXDw8PHAAx39wrGjh0BScbkbCQm6zljXBKTdZ5ETCtWrDDmzZtnzJ492/jhhx+SXd5y3LTqGCtXBqlkGLBmjQeTJkHu3LBxI3z4oa33YbBo0SL69+9PTEwsLi7+zJzZierVbTxCTQiRLjRr1sxu25ZxBqlw6RLUqwfDhxehQQOzrpCtE0FwcDBvv/02Pj4+1KlTh3/+OcagQV3p0EESgRDC9iQZpNDatVCmDOzfDzNnnmPJEsif33bbv3XrFoMHD6ZcuXJcunSJzZs3s3z5cooUeZaxY9PeK0kIIR5FkoGVIiKgXTto1AiqVjWvBmrXTtmk2PcYBgwZYn5PaOPGjZQqVYqpU6cydOhQAgMDHTYARQiRucgzAyvs2QOtWpn9l7/5xuyz7+ICV6+mbnuJxwycO3eOXr16sW7dOmrVqsWmTZtQStn2lxBCAOZJ2OjRMHx46tZv06YNAwcOJEuWLKxYsYLr168TFxdH48aNWbJkCbdv3+b27dsMGDDA7pPs2JJcGSQhLg6+/BKqVzfrCB05Yh6803Kr5sExA7GMGTMFT09P9uzZw3fffceWLVskEQhhR/Pmmclg1arUrZ81a1ZWrFgBQIECBejZsyf9+/dn8+bNAIwePZoRI0awcuVKW4X8REgyeIxz58wyEn5+5lSUu3bZpu/y/8YM7EFrL4YP78+nn37KsWPHaNmypVV9l4UQqRMRAdOmwZ07MH48xMSkfBs5cuQgKiqK69ev8++//1K8eHGKFStGvXr1ALMX4I0bN8iRI4eNo7cvSQaP8PPPUK4c/Psv/PGHWQ/IMgAwzRo3voqHhy9QFXf3rGzdug9/f/8nPvm1EJnRrFn3Tsbg0KHUV/Jt3Lgxa9asId5S8yI2NpYLFy4A4Ofnx6JFi2jTJn1N1SLJIIHYWHP2sXr1oFIlOHwY3nnHNtuOj49n4cKFvPmmIiZmJS4uM5k5cz81arxpmx0IIZLVvj289JL52tMzdTW7XFxceOutt9BaU7hwYS5cuMC///7LkSNHAPjiiy8YP348Tz/9tA0jtz9JBhbnz8Pbb8PUqebD3Z9+Alv9LQMDA6levTrt27fn3Xff5Z9/NIMG9ZAxA0I8YYULQ/fuZiG6vn3NAaOp1axZM27cuMH48eOZNm1amievdzTpTQRs2wbNmkG2bLBjh9l11BYSTkRfvHhxtmzZcv8fjBSPE8Ixevc2T/5SW8l31qxZALRu3fqB90NDQxk3blxaw3OYTH1lYBgwaRLUrg2lSsHBg7ZJBIZhsGbNGjw9PfH398fPz4+jR4+m+zMHITICV1fzDoD01XhQpk0GN2+aVwMDBpiXi1u22GYayn/++YePPvqIxo0bU6ZMGYKDgxk2bBjZbT29mRBC2FCmvE109qw5If3x42Zf40TzUaRKdHQ0EydOZOzYsRQsWDDZeQaEEMKZZLpksHs3fPKJOQ/x7t1mF9K0+u233+jevTunT5+mT58+jBgxgtxpeTIlhBBPWKa6TbRkiTmQ7JVX4O+/054IwsLCaNKkCe+//z5FihTh8OHDTJw4URKBEE7IiDcIPxP+0JcRbyS/soNduHCBsXbudZIprgzu1SIZMQLatjVHAaflFn5sbCwzZsxgxIgR5MmThyVLltC6dWu5JSSEE4s4F8GMF2c89H6v073I94LzD/q09/ElwyeD2Fjo0sWcmH7kSPMrLZ/pjh076Nq1K6GhoTRv3pzZs2fL6GEhMpEmTZrw1FNPUahQIU6fPs1zzz1HZGQkTZo0IWvWrEyePJksWbIwaNAgYmNjmThxItmyZaNHjx5s2rSJ2rVr89xzzzF27Fi6d+/OuHHjyJ49O507dyZHjhzMnj2bqKgofH19cXd3Z/LkyXh4eJDflrXyHyFDJ4OoKPPh8JYtZjJIbb9iMG8J9evXj++++46KFSvy999/kzNnTkkEQmQyefPm5f333ydr1qxcvHiR/PnzkyVLFk6dOkXJkiVp2LAhW7du5eDBgxQqVIjcuXPTtGlTihQp8sjthYWF8cMPP+Di4sLEiRPJnj07OXPmJDAwkDNnzjB06FDCw8PZuHGjXX+vDJsMwsPN2ceOHjWnpHzvvdRt5+7du/j7+zNixAiyZMnCvHnzaN++Pa6uroSGhto2aCGE03N3dwcge/bsGIZB165duXLlCocPH2blypUUKFCAqlWrEh8fzwsvvICPjw8bNmzg7NmzuLi4EBsby/nz5+9vz8vL6/4tIBcXF1q3bk3+/PkJCgrizJkzZMmShay2Ko6WBJsnA6VURWCipdkWmAtEAr6AO7D4cW2tdaQtYrhyBd59F06fhq1bwds7ddvZvXs3Xbt2JTAwkA4dOjB27Nh0V29ECGFbCe/du7m5MWjQIOLj42nRogXx8fEEBAQQHh6Ot7c3169fZ8GCBURHR1OrVi2KFCnCsmXLKFq06EPbAmjevDkTJ04kPj6eHj168MknnzB16lSyZ89u944pLkbi6bbSSCnVClgDvAssAqphHvQrAy9gHvyzP6qttZ71uO0GBAQYXl5eye7/wgWoUwf+7//M20Ovv57y3+G///5j4MCBLF68mPLlyzN79my8H5FRQkND8fT0TPkO7Ehiso4zxgTOGVdGicmIN4g4F/HQ+x7Pe+DimvaHs874OQUEBODl5WXVL2fzrqVa6++01tFABSBAax0MBAIvAYW11kFJtNPkwgVzIprISPjzz5Qngrt37/LVV1+hlOLHH3/E39+f/fv3PzIRCCHSFxdXF/K9kO+hL1skgozALs8MlFKvAkWB05a33IDcQHQS7VzJbTepe/TXrrnRtu0L3LzpypIlZzCMWFJyS//gwYOMGjUKrTWNGjWiT58+FChQgOPHjz92nejoaKd7biAxWccZYwLnjEtiso4zxpQS9nhmkAPzmcFnQDullBfmbaCzQFQS7XPJbftxl2Dh4eYcxTdumFcEr732itXxhoWFMXDgQJYsWcIbb7zBX3/9RaVKlaxa1xkvCyUm6zhjTOCccUlM1nHGmAICAqxe1h5XBoOBgpgJYSfQG/MqoCMQA/gn0U6xmzfho4/MWcm2b4fXXrNuvdjYWPz9/Rk5ciRubm58/fXX+Pr64uYmcwwIITIfmycDrfUIYESCt5YmWiTxXHCpnhsuLs4cR3DkiPmw2NryEtu2baNHjx6EhIRILyEhRLrUrVs3/P39bba9dF2baNgw2LQJfvgBKldOfvnz58/TvHlzatasSe7cudm/fz/z5s2TRCBEZmIYMHy4o6NIM1uXp0i3g85WrYLx42HcOPjgg6SXjYmJYcqUKYwZM4ZcuXKxYMEC2rVrh6trus6FQojUWLcOpkyBmjXNypUplJZyFE899RRDhgwhZ86cNG7cmLCwMPLmzUutWrXun+m3adOGokWL0rBhQzw9PZk4cSK3b9+mbt26lC5dmuHDh/Pss8/a/GNJl0fDo0fhs8/MW0QDBya97M8//0ypUqUYMWIEHTp04Pjx4/j4+EgiECIzMgz4+mu4fRtmzkzVJu6Vo6hYsSLZsmUjf/78FChQgFOnTuHh4UHDhg3JkSMHBw8e5PLly+TOnZvWrVtTpEgRbty4QXx8PI0aNaJkyZIPbPfemf758+cZOXIk3t7e/PTTT0RHR1O4cGGCgoL49ddf+fTTT+nevXuaP4rE0t0R8do1+PhjKF4cFi58fNE5rTUffvgh9erVo1ixYhw+fJgZM2ZILSEhMrN16+CPP8zXv/5qToCeQo8qR+Hj40O+fPlYuXIl165de6gcxR9//MGmTZvIkycPvXv35tixYyxevPh+eYqIiAiioqIAKFmy5P19uLi4UK9ePXr37k3x4sUB7pensPWA4XR3m2jYMLh6FTZvhkeNzo6MjGTUqFHMmDGDZ599llWrVtG4cWMpLy2EMB8wJuxyuHRpim8VpaUcxZ07d/jmm28wDIM6depQokQJpk+fzpkzZx5IAPd89NFHjBgxgjVr1tC0aVO8vLyYMGECBQsWtP0xzTCMdPF14MAB4/x5w8iWzTBGjTIeEhcXZyxcuNAoXLiw4e7ubowcOdK4efPmwwvaUEhIiF23nxoSk3WcMSbDcM64JCbrOGNMBw4cMAwrj7Hp6jbRxImQIwf06PHg+3v27KFixYr4+PhQrVo1QkND8fPzI2fOnI4JVAgh0pl0lQzmzYOePcHDw2yfPXuWli1bUrVqVeLi4tixYwerVq3ixRdfdGicQgiR3qSrZJAlC/TuDVFRUQwbNgylFFu3bmXevHkcOHCA6tWrOzpEIYRIl9LVA+SuXeNZs2YBI0aM4Nq1a/Tp04chQ4aQN29eR4cmhBDpWrpKBmvXvsXJk3/RtGlTxo8fz0svpbnqtRBCCNJZMnj+eXe+/34/b775pqNDEUKIDCVdJYOtW7fKeAEhhLADm097aS8BAQHpI1AhhHAi1k57mW6SgRBCCPtJV11LhRBC2IckAyGEEJIMhBBCSDIQQghBOuhaqpR6BlgMRAK+WutIB8fTBsgHrHSGuJRSFYGJlmZbYK4TxOQBzACyAXOAwY6OKUFsPwCTgC+cISallDcwEjgNLHKiuEoDrYGpOMe/8xZAdcwT2OcBFyeIyRvoBeQCeuAc//fqWGJxAdqTgr+d0ycDoL/lKzvmwW6WowJRSpUF3gaOOlFcrwLvA+8Ch4FqgLuDY3oTGA54YH5WZZwgJpRStS2x9cM5/nYA5YD+WutgpdQkJ4rLD/NgMswZYtJaLweWK6XqA5OBj3H8v6lWQDegLHAIM1k5OqbOWuv6SqmPgEuW2KyKKT3cJiqstQ4CAgGH1p/QWh8FlmBmXaeIS2v9ndY6GqgABGitg50gpt8x/20NAtY5Q0xKqSxAO8wzpRhn+NtZvAH0UUqtB8o5Q1xKqcZASeBb4GVniCmBhsDfzvBvCvgN8wDbE9jkJDG5KKVyWGKIT0lM6SEZRFu+u2FejjnavQEcThOXUupVoCjmrQZwgpiAq0Aw5pkJOD6m7sBMIB64N5u4o2MC6Ke17gCMA96xvOfouFoCHwC+QH7Le46OCaWUAk4CtyxvOTqmSpi3i1cDlS3vOTqmScB8YCjmiStYGVN6uE0UpJTywrxMPefoYBJwirgsZwETgc+Adk4SUyNgO7Ac8FNKlce8VHXk368cUAzzP20FJ4kJzAPvPKAAMMIZ/n7AdeCCJQ4vJ/qsPsQ8G3/LST6n0sAoIC8w2xk+J631X0qpXMBBIC4ln1N6SAbfAv5AbqCjg2O5x8CMazaOj2swUBAzIewEejtBTMcwz0puY96+GuDomLTW7QCUUiMscTg8Jotoyy2iu0AfYCyOj2sOsBTIA9S1xOXomAC8tdbTlFIncI5jwiLMk554oD7O8zm1Bz7FTExWf05SjkIIIUS6eGYghBDCziQZCCGEkGQghBBCkoEQQggkGQiRbiil3CxlGVKz7htKqTK2jklkHNKbSDgNpdRIoCoQBpzUWn+plHLDHNRzFZiptT6qlOqGWRZkm9b661Tu6zmgl9Z6QBriHYBZguAu0FRrfSe127Jyf4OAVZilPiIx6xj9H/AisB/4HbNuVi+t9T+J1nUBpmOWvrBrnCJ9kisD4WxGa63baq2/tLSLAOFa6w6WciAAn2itm6YkESiltiVsa60vpCURWDTQWtfFHFPx2qP2YyuWwYVuWut/LW+9CryOOZYDIERr3Rqz6F3XxOtrrQ3McTHN7RGfSP/Sw6AzkbkMVkq1wzzL/RtzII+XUmosZtG0QUBxS3sdEIpZWdMALmuthyulxmGWnLhXj2gU8JJlnYVa65NKqRcAP631Z5YRmwu01s2VUn2BAMwrkT6Yw/i/01pveESse5VSi4BRWut/lFJDEuxnndb6b6XUdmAvEKm1HquUao5ZdqIw8DnmCdkkIAr4WWu9Uim1ALPcSQ6ttY9lX56YxdDuMTBLfRzGLJFyr0yKBxCplHoFs3JsFLBWa71Sa62VUs2s/kuITEWSgXA247TWf95rKKUmAO5a6yGWt8YqpWrfayulOmOWcjgOlLHUr8mntW6nlCqHeUAdrJTyTrCNewwArfVNpdRNpVQ+oIrWeopSagVwETMpVAEeSgZa675KqZLAXKXUQMvBvnai/TwLLNJaH7O0hwA/YyaBNzErS8YBCzFrOQEozAN5UKJdJpzY/CxQAzhhab+mlPoWuIM5ujo7Zknl8kBjzBo6939nIRKT20TC2bhY8V7i9nrLAfinRD97DvMM3Jp9/YB5ED1uaRuYVxzDMas+PkAplV8ptURrHYJZubKa5UdxSqmE/68uJkgEADcssU7FrP9zCvgSMwGMtTwj6Q2EAwssVy1gJgqvBNs5AhTCLIUAEKq1/kxr3UlrHYFZZO51zOctWSwxl8Ys9CbEQyQZCGczRCm1RCk1L8F7ic9mE7a/A+oppZYCz1oOvDeUUkuAFsAZy3K7lFLLEvXGSbidLZj1XL6ztMcBX2EW27uWOEit9XXgP8sVRDvLcmAefH9QSvV8TOz+SqlVwATM4mEFMOtLVQUOaa3jgM7AJ5hXJrct+4vBvP1TyrLNI5hzRTzOGcwE0x3IZSnh3Qkz6QnxEOlNJEQ6opSq/5jnF8mtVwa4lbiXkRD3SDIQQgght4mEEEJIMhBCCIEkAyGEEEgyEEIIgSQDIYQQSDIQQggB/D+3NCcfQmFQDgAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x10bea160>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "fig_bowers, ax_bowers = plt.subplots()\n",
    "# draw virgin(loading) curve\n",
    "ppp.plot_bowers_vrigin(\n",
    "    ax=ax_bowers,\n",
    "    well=well_cug1, \n",
    "    a=a,\n",
    "    b=b,\n",
    "    vel_log=vel_log_filter_smooth, \n",
    "    obp_log='Overburden_Pressure', \n",
    "    upper='T12', \n",
    "    lower='T20',\n",
    "    pres_log='loading', \n",
    "    mode='both')\n",
    "\n",
    "# draw unloading curve\n",
    "ppp.plot_bowers_unloading(\n",
    "    ax=ax_bowers,\n",
    "    a=a,\n",
    "    b=b,\n",
    "    vmax=4600,\n",
    "    u=u,\n",
    "    well=well_cug1,         \n",
    "    vel_log=vel_log_filter_smooth, \n",
    "    obp_log='Overburden_Pressure', \n",
    "    pres_log='unloading')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 3. Pressure Prediction with Bowers model"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "predict pressure with coefficients calculated above:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "pres_log = well_cug1.bowers(\n",
    "    vel_log=vel_log_filter_smooth, a=a, b=b, u=u)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "View Bowers Pressure Results:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAARUAAAHyCAYAAADFr1rYAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAIABJREFUeJzs3Xd4FOXexvHv7KYn9AChJZTAEDpp9CZ4PMgREbCiKIqiiBSVI3oErMfee0VQEI+votgbXYFA6BAGkB5Ch5CQnuz7xyT0JJtkdmdn9/e5rlwyu5uZGyQ3s7PPPI/icDgQQgij2MwOIITwLlIqQghDSakIIQwlpSKEMJSUihDCUFIqQghDSakIQ6mqulBV1dpm56gqVVVfU1V1qNk5rEhKRRjN8gOfVFWdCvQ0O4dV+ZkdQJRPVdXpQC/0/18pwLPA28AeIEDTtNGqqrYG3gGCgWmapv2qqmpv4Dn0fzwe1zTtR1VVrwPuL971vcBR4E1N065SVbUPMEzTtPGqqm4EXgfGaprWWVXVEOBToC7wq6ZpT5USVwHeUFW1JfC6pmmfqaqaCLxc/NyjQCQQCPQA5gN3AkOAO4CbgAzgOk3TTqmquhCYBYwHrtA07bCqqnOAhsBxYLimaUWX+DOrB8wEagI/apr2ZPG+hmmadlxV1Y2aprVXVfVGYDJwDMgErgU+Ls4qKkHOVKwjSdO0vkB9oAbQF/hC07TRxc8/CzwIDAAeK35sOPA4cBlQUPzYCOA24BrAv/ixS51dBAHhmqZ1Lt6+C/hL07TeQP9y3uJ8hl6CDxdvvwoMBQYCLwFrABU4DLQG7EA1YJSmaT2A/yvOWCJO07TOxYVSG2hW/GfxFhBWSoaH0cuyGxCmquqFf9dLfs/Xof+Z/A28o2lagaZpqUipVJqUinVsKP5vCvqp+UpN0xaf83w08ALwHRBS/NgbwO3Al5wtlafRC+hdIOuCYwSe8+sC4MVztlsBNxb/a18b/UyhNGs1TcsFTqiqGgqEapp2WNO0U0AOsLt4f7lAY/RyaQY0UFV1ATCq+Bglni/5haZpx4GvVFX9Gb2oLvw9lFDRywtgyiXOZkp+r2vRz/DqAkvL+D0JJ0mpWEfJGUM7YCFw+oLntwE3aZrWD/30HaA/+g/oncD04sd6a5o2BP0H9QEgn7M/wH3P2V+Bpmn552xrwAvF+38a/e1CaTqqqhoI1NI07TSQoapqfVVVa6C/XctA/yHejV4m64FdwBpN0y4rzrzynP1llvxCVdXGwH5N0/4JBKCfEV3KViC2+Nd/qKpaH8gDahW/Vaxf/Fy8pml9NU0bpmladhm/J+EkKRXr6KCq6iJgL/q/zhe+ZfkPMFtV1RWcfVuzD/gD/brFJ8WPnSx+zevA55qmHQDSVVV9C/1tVYkL9/8hcK2qqkuBwehnF5fiQD87WsrZM4z7ga+BH9GLDPQzrx3oxbJO07QjwC+qqi4BZgD7S8lxCBiuquqfQAtgdSk5ngPGF+/vV03TDgFz0c/ebgEOFL/uiKqqq1RV/UlV1TFl/P6FkxS5S9nzqao6DdikadrXZmcpUfwv/1z0Hz6l+L8/aJr2kglZFnJ+CSRrmjbZie9T0As3AL2ot2ua9m/XpPQdUipCCEOZ/pHyOR/9nQLuLL6YJ4SwKNNLBX2MwGT0q/EjgTfNjSOEqApPuFAboWnaJmAj+icBQggL84QzlZzi/9qB0NJelJycLBd/hDBJXFyc04MBPaFUNqmqGof+9mdfWS+Mi4tzTyIDpKSkEBMTY3aMCrFaZqvlBWtmTk5OrtDrPaFUZnB2uPVdJmcRQlSR6aVS/GnPLWbnEEIYwxMu1AohvIiUihDCUJYqlSLHRdNmCCE8jKVKZdPhTWZHEMIpDz/8MJmZmaSmpvLf//73kq9JSkpi5syZbsmzYcMGvvnmG7ccy1KlsmDXArMjCFFhilL6EI+ynitLSWk5+3yHDh0YMmRIpY5VUaZ/+lMRR04fMTuCEE577LHHyMvLIyIigkOHDnHgwAF+/fVX0tPTmTVrFjabjcTERJKSkvj666/Jyclh2rRpHD58mLfffhu73c6UKVPIysri+eefJyAggPvuu4+NGzeyadMmXn31VYYPH054eDiPPvooAQEB3HHHHezcufO851u3bk1SUhIpKSnceuutZGZmMn36dAoLCxk4cCBXXHGFob9vS5VKaEAoJ7JPUCu4ltlRhBdYtHsRi3Yvuujxvk370rdp33JfX9rrSjz22GOkp6czbdo0fvrpJ44cOUKvXr2YO3cur776Kr///jtHjx4FIDMzkzff1G97e/rpp3n++edJS0vj008/pVevXoSFhXHdddfRqFEjoqOjSUpKYuLEiYSFhZGRkcHw4cNJSkrir7/+4p577jnv+RIlZ0Xz589n6NCh9OjRo8ID25xhqVIZ0noI7yW/x5SeU8yOIrxAeaVQ1deXaNGiBRs3biQoKIhatWqhKAo2m42AgIAzr0lMTDzza0VRKCoqomRakqioKG6//Xbmz5/Pvn37GDp0KHa7ncLCQgD++OMP9uzZQ58+fVi7di3Aec9f6NzpTlJSUgwfqW6payrZ+dmsSl1FXmGe2VGEqJAmTZqc+eG99tprmTp1KvPnzwcuvq5y++23M2XKFF577TVGjBjBiRMn+Oijj9i3bx+RkZEA9O7dm0cffZSvvvqKevXqsXfvXubPn09GRsZFz19o8ODBfPXVV0yYMIHatV2wRJPD4bDE1+rVqx0FhQWOu+bf5Zi1bpbD023ZssXsCBVmtcxWyZuUlOQYN26cIzs72zKZz7V69WqHowI/q5Y6U7Hb7NQNrcu6g+vOO4UTwpMlJCTwxhtvEBQUZHYUt7BUqQDUDKpJl8Zd+H3n72ZHEUJcguVKpWdkT4L9gvnl71/MjiKEuATLlUpcgziKHEW0qdtGRtgK4YEsVyr+dn+ubn01N7W/iTkb55gdR4hLunBE62+//cby5csv+do333yTrVu3ujTPoUOHeO+991x6jBKWK5USQX5B1AyqycHMg2ZHEaJcl19+Od26dTuzXd4w+/KUV0QXPl+/fn3GjBlT6uuNZKnBbxe6vfPtvJ/8Po/0esTsKEKUad68eVSvXp2GDRsyb96884bRA3z55ZdkZmbSo0cPEhISePbZZ1EUhbvvvpv69evzyCOPEBISwvDhw8nOzmbZsmWkpqZyzTXX0LFjRx588EH8/f258sorsdls5z2fmJhIamoqM2fO5JFHHqGwsJCpU6eSn59Phw4duOUWY+dIs3SphIeEk5Wfxem804QGlDpnthCX9M3Wb1h3cF2lv79TRCeGtK74TXrXXHPNecPof//9dxITE0lISOCFF15g27ZtTJw4kYiICB599FEmTJhAUVERw4YNo23btlSvXp0tW7YwYMAAWrduTW5uLsOGDSMlJYVFixbx1FNPnfd8iZJBdsuWLaNTp05cd911rF5d2qqxlWfpUgEY2XEkn274lLvj7zY7irCYIa2HVKoUnFXWWKoLh9FHRUXh7+9/ZtvhcFBUVISiKFSrVo2JEyeyfPlykpOTue+++7DZbBQV6fMLrV+/ngULFjBixAg++eQTgPOeL0tKSgrx8fGV/B1emmVLZemepcQ1jKNVnVZ8kPwBhUWF2G12s2MJccb06dPx8/OjW7duFw3FLxlG37dv34u+b8SIETzzzDNn3v7k5eXx4Ycf4nA4uPzyywHo1q0br7/+OvHx8fTu3ZuTJ08yZ84csrKyLnp+9OjR5+2/R48eTJ06laSkJNq3b2/479syayknJyc7zr3xaeGuhRQ5iujfvD9L9yzlaNZRrom5xsSE57PiUgxWy2y1vGDNzMnJyRVa98eyn/50adyFFftXAPqAuKV7l8rQfSE8gGVLJcQ/hKx8/VRPURT6RPVhyZ4lJqcSQli2VAAC/QLJKdBXTb1KvYrvtn1nciIhhKVLpUujLiSlJgFgU2y0q9eODYc2mJxKCN9m2U9/AHpF9eJY1rEz2ze1v4lpC6fRoX4HE1MJoY+Yzc3NJSMjg3HjxtGxY0ezI7mNpUslxD+EkBohZ7YD7AE0rNaQ3Sd307RmU/OCCQE88cQT7Nu3j99++01KxcpGdRrFi3+9yOP9Hjc7ivBx06dPZ9OmTTz77LOMHz8eu93O0KFDeeihh2jevDkNGjTgzz//5JZbbmHOnDmcOnWKhx9+mNTU1PNm109PTz9vNv3o6Gizf2tl8rpSqRZYjUC/QI5mHSU8JNzsOMKDffMNrKv8KH06dYKyltJ5/PHHOXHiBE8++SSvv/46aWlpvP/++4SGhuLn58eePXsIDw/n448/pkGDBvj5+ZGSkkL16tXPm11/+/btZ2bTb9iwYeUDu4nXlQrAHZ3v4KM1H/FQz4fMjiI82JAhZZdCVTkcDkJDQ8nNzT1vdvyIiAjS09NJTU2lX79+pKWlMXbsWI4cOUJ6ejq5ubnnza4fGRl50Wz6nszypZKfD35+cO4o6Pph9cnMy5QbDYWppk6dSlZWFpMmTWLKlCkoisKwYcNQFIXNmzcTHBxMmzZtiIqKYsqUKRQVFTFt2jT27Nlz3n5OnjzJRx99RE5ODv379zfpd1MBFZkl28yv4hm9L/L99w7HhA8/c+QX5p/3+I5jOxxvrnzTuenCXcCKs6ZbLbPV8joc1szs1bPpX8qxY9Aiot5Ft7C3qN2Cvel7yS/MNymZEL7J8qWydSvc2K8zS/csvei569tdzxebvzAhlRC+y/Kl4uenT9Z0PPv4Rc/FNohlbdpaudFQCDeyfKmUXKBVFIUix8WT0gxsOZCfdvzk5lRC+C6vKZX29dqjHdUuer5/s/6y8JgQbmT5UikxvM1wYupePPmNoij0aNKDZXuXmZBKCN9j+VIpmYbzwun6znVNzDXMS5nnpkRC+DbLl0pgIJS3fIpNsdEpohNr0ta4J5QQPszypXLFFfDrr+W/7oZ2NzB301zXBxLCx1m+VJo3h927y3+dv92f6NrRpBxJcXkmIXyZ5Utl1SpISDi7vefknlJfO7LjSGatn+WGVEL4LsuXSmoqREae3f5gzQelDnYL8guiQbUG7Dqxy03phPA9li+VunXh8OGz2y1rt2TH8R2lvv72zrfz8dqP3ZBMCN9k+VKpV+/8UunTtA+L9ywu9fVhAWFUC6xGWkaaG9IJ4XssXyr790OjRme3o2pEsfvk7jK/587YO/lgzQeuDSaEj7J8qezYAS1bnt1WFOXMvA6lqRVcC5tiO28mfiGEMSxfKllZEHrB5G59m/Y9s3phaeRsRQjXsHypXGp0/uUtLi93Gsn6YfXJys/iVO4pFyUTwjdZvlSqYnTsaD5a85HZMYTwKpYvlarMvxRZI5IjWUfIzs82LpAQPs7ypRIcrF9XqazbOt3GJ+s+MSyPEL7O8qUSHa1/AlRZreq0YtfJXeQV5hkXSggfZvlSadUKtIsnfAPg+23fO7WPmzvczGcbPjMwlRC+y/KlEhMDq1df+rmV+1c6tY8O9Tuw+fBmWc5DCANYvlT8/aFFC9i48eLn6obW5VDmIaf2M6LDCOZsnGNwOiF8j+VLBWDUKJgx4+LHe0f1ZsmeJU7tI7ZBLBsObaCgqMDgdEL4Fq8oFX9/6NEDFl9wH2H7eu3ZePgSpzCluKHdDXyxSRYfE6IqvKJUAK65Br6/4Lqs3Wa/5FpApUlolMCatDUUFhUanE4I3+E1pWKzQfv28Oef5z9+e+fbK7Sf4W2G839b/s/AZEL4Fq8pFYBbboH//Q9ycs4+1rxW8wrto1uTbqxMXVmhMxwhxFleVSqKAhMnwvvvV20/Q2OGyjpBQlSSV5UKQLNmcPRo1fbRM7Iny/Yuk7MVISrB60oFLj0dQkVd3fpqvtO+q/qOhPAxXlkqRugT1YdFuxeVOYOcEOJiXlkqF56pFBYV8tyy5yq4D4VBrQbx4/YfDUwmhPfzylK5kN1mJyMvo8Lf179Zf37f+bucrQhRAV5ZKpfqgLCAMDJyK1YsiqJwRfQV/Pq3E4s1CyEAHyqV7k26s3z/8grv64oWV/Dzjp/lbEUIJ3llqRReYpR9YqNEp6dCOJeiKPRv3p8FuxYYkEwI7+eVpdKgARw4cP5jQX5B+Nv9K7W/QS0HOT3hkxC+zitLJSrq4lIBmNJzSqX2pygKfZrqHzELIcrmlaWSl6dPh2Ckq9Wr+Xbrt3JtRYhyeGWpJCTAL78Yu09FUejbtC8Ldy80dsdCeBmvLJXISH0p1ORkY/c7WB3MfG2+nK0IUQavLBWAe+6BTz6BfAPnslYUhQHNB/DHrj+M26kQXsZrS8Vmg7vvhg8/PP/x03mnOZh5sNL7HdRyED9s+0HOVoQohdeWCkDbtpCaCidPnn2syFHEzHUzK71PRVG4vMXl/LbzNwMSCuF9vLpUAO67D1577ex2tcBqlboP6FwDowfy4/Yf5WxFiEvw+lKpX1//iPlc/jb/Ki0cpigKA6MH8vOOn6uYTgjv4/Wlkpt7canENohlTdqaKu33Hy3+wa9//ypnK0JcwM9VO1ZV9RagJvAFMBM4BdwJBJW1rWnaKSNzzJgBd9xx/mM9I3vy0dqP6NK4S6X3qygKV7a8kh+3/8igVoOqmFII7+GSUlFVtQPQB9gATC7+CgRGAlHlbL9pZJbUVGjd+vzHagXXIiY8psr7HtB8AJN+mcSVLa9EMWIOSyG8gEve/miatgGYBShAhKZpm4CNQDMntg1lt1/6cSPOLhRFYbA6mO+2yVy2QpRw2dsf9EIBKFmFxw6ElbMdWtYOU1JSKhziyJFwUlKqOL1+GSIcETy//nmiC6PPO1vJycmpVF4zWS2z1fKCNTNXlCtLpcQmVVXj0N/e7AUyy9jeV9aOYmIq/palbl2Iialb4e+riNtCbmNbzjaGtB5y5rGUlJRK5TWT1TJbLS9YM3NyBe93cXWpOIAZwNvoZyV3AbnAW2VsGyo31+g9Xqxv075M+nkSg9XB2BSv/0BNVEFIwB5OH9tMUdFp/AIaEVS9H1nH5xFU4zLsfrUpKjxN9smfUGyBBIZ1x+5fx+zIFeayUtE0bTGwuHjzlgueLm/bMK1aQVISJCa66gi6a2KuYV7KPIa1GebaAwlLy8qLIjAsjILc3QRV70VOxjIcnF20LjfzL4Kq98LmV4fs9D8IqXmFiWkrx+v/Wb3tNpg9G9LTL35u0+FNbDu2zZDj9I7qLasaigoLqtYTm736me2iwkzs/vVRFD8UxR1XJ4zn9aWiKDBtGkyfDkUX/LyHh4QbOlP+8DbD+WrLV4btT/iec4vEUZRTxis9l9eXCkCdOnDrrfDyy+c/HhEWUaU7li/UI7IHy/cvl7MVUWk2v7oU5KXiKMrFwSVmcLcAnygVgM6d9Qmxf7zEgoNGDrW/ts21fLn5S8P2J3xLYGgn8rLWk3XyJ4LCKj/i20yWetOWk/FXlb5/2GB47IkmtGp+iMaN9BuCWoTaWb/nc1rXaWpAQuhcE+au/4bWkZFVzutuRYX+lspstbygZy7IOwGK/Ux2m706+dlbKbnF1WavCXbIz/mb/Jy/zQt7RmCFXm2pUgmq1r3K+5j+BDz2WBOeK15aeUDrJvyw/Qc6Na36vkuMjA3kp/W/0bFD5WbvN4vNnkJQNeuMobBaXrBmZqjYOBWfeftTIjgYAs8p3iY1mjCopbE3BMY1jEM7qZFXmFf+i4XwMj5XKnDxp0BNajQx/BjXNLuGT9Z9Yvh+hfB0PlcqDoc+f62rtazRkr+P/012frbrDyaEB/G5UlEUY2fYL8udcXfy0dqP3HMwITyEz5UKQEQE7Nzp+uNE144mLSONzLxM1x9MCA/hk6Vy113w9tvnP+aqaSHHxI/hvdXvuWTfQnginyyV4GD9BsMlS84+NnP9TPam7zX8WJE1IknPTedkzsnyXyyEF/DJUgG49lr46isoLB4J3a1xNxbvXlz2N1XSPfH38M6qd1yybyE8jc+WiqLAyJEws3hdsVZ1WqEd01xyrAbVGlBQVMCR00dcsn8hPInPlgpAXBxs367/2tUTV9+TcA/vrn7XpccQwhP4dKkAtGhxtlgahDXgQMYBlxwnPCQcf7u/y/YvhKfw+VIZPhy+LL6puG/Tvqw/uN5lx7onXs5WhPfz+VKpWRNOndKH7ret15aBLQe67Fg1gmpQM6gmu0/udtkxhDCbz5cKQPv2sM2YWSXLNSZuDO8nv++egwlhAikV9HEr7ph1HyA0IJSG1RoaNjeuEJ5GSgUICYGTbhybNjp2NB+u+dB9BxTCjaRUgP799Yu1eW6a/iTIL4iWtVuy8dBG9xxQCDeSUgH8/WHSJM7MBrf16FaXH/PWTrcyc/1Mlx9HCHeTUinWooX+tWgR/Lj9R45nH3fp8QLsAXSs35FVqatcehwh3E1K5Rw33gjz5kFi/V4s3bPU5ce7qf1NzNk4x+XHEcKdpFTOoShw332wYHYsaw+udfnx7DY7vaN6s2DXApcfSwh3kVK5QHQ0dO5kZ81P7d1yvCGth/Dt1m9dNp+LEO4mpXIJV10Ffo5Qflvk+hnbFEXh6tZX8632rcuPJYQ7SKmUYuLI5qxMds9nzJc1u4xFuxdRWGTNZS6FOJeUSmmOt6JLh9puO9yI9iOYvXG2244nhKtIqZRi3z6IjHTf8RIaJbDh0AZyC9x0v4AQLiKlUorGjfVicafRsaNlSQ9heVIqpUhIgBUr3HvM1uGt2X9qPxm5Ge49sBAGklIpRUiIPs9KoZuvnd4df7dM5CQsTUqlDDfcAPc8udqtx4ysEUlWfhbHso659bhCGEVKpQyxsbDnYDorV+e49bhjE8by1qq33HpMIYwipVKO+x8o4sV3jpKW5r5j1g2tS6A9kP2n9rvvoEIYREqlHD2jutPuhs956inIcOP107EJY2UBMi9UsHMx2z9/kC0fjGLvTy9TVJjP9s8nk3NM/6jx1K7VbPt0Ats+m8ipXe59620UKZVyhAaEkm9P57HH4JFHID/fPcetFliNBtUayLSTXsaveR8iut9MjZbdiRx4PweXzsRRePYvVe6xfbS86SWaXvUwR9fMNzFp5UmpOMHf5k+N2nnccAPMneu+48q0k96vYd/RBNSMOLNdN/4aFLsfGXvXE9akg4nJKk9KxQn3Jt7Lxg0Kc+fC1Ve777hBfkHEhMewJm2N+w4qTJd36jAntiygbvwQs6NUipRKOQoL4bP3w/lzqT+vvgrVq7v3+Ld0vIVP13/q3oMK0xTkZLD3xxeJHHg/is3P7DiVIqVSjmefhQEDYPx4sNvdf3w/mx89I3uyaPci9x9cuIXC2XW89//+FoV5WaQtm8XBPz8zMVXlKVaZHCg5OdnR4PQGtx93zvctad38BLFtjlbo+44GRROes8OQDA6Hg+k7Pufx6BtdupC8kZndwWp5wZqZ00I7EBcX5/RfPEudXzXsPcrtx3ygl/6pT8t+oKrOf196SgoNE3sZluOGJpEk52cxWB1s2D4vZHRmV7NaXrBm5rTk5Aq9Xt7+lENR4Mkn4Z13YH9qkWk5BjQfwIJdC2QiJ+HxpFSc4OcHzzwDQ+9dxRoTP4i5ucPNfLbBmu+zhe+QUnFScDBcMe4n1q4tYsoU2L3b/RniG8az6fAmsvOz3X9wIZwkpVIBnRq2J/bK9Uybpq8P9MYb7s8wJn4M7ye/7/4DC+EkKZUK6BXViyV7lhASAmPGwObN7s8QXTuaY9nHOJF9wv0HF8IJUioVUC+0HodPH2bhQv0ToalTzckxLnGcTI0gPJaUSgUcOwYrZw7h8GF45RVo1MicHPVC6+Fv82dfupsn0RXCCVIqTnI44D//gf+9msD11+sfNZtpbMJY3lktUyMIzyOl4qSjR6FdO6jtvqWAylQtsBqRNSLZfNiECztClEFKxUmHD0O9emanON8dne/g47Ufmx1DiPNIqTjp6FEIDzc7xfn87f50adyFpXuWmh1FiDOkVJx05IjnlQrA8DbD+b8t/4dVbgwV3k9KxUmnT0O1avqvC4sK2XJki7mBitkUG0NaD+Fb7VuzowgBSKk47fRpfYExALvNzucbPzc30Dn6NevHkj1LKCgqMDuKEFIqzjp48PwLtQ4cHvWW49aOtzJz3UyzYwghpeKsgoLzx6Y0r9WcXSd3mRfoAh0jOrL9+Hay8rPMjiJ8nJSKk9q2hV9+ObvdO6o3i3cvNi/QJcg6zMITSKk4acQI2LABlizRt1vUasHfJ/42N9QFmtZsSkZuhqzDLEwlpVIBDz4If/4Jf/0FiqLQv1l/syNd5N7Ee3kz6U2zYwgfJqVSAYoCU6bAr7/Cjh36py6eJjwknLCAMHaf3G12FOGjpFQqSFFg4kT45huzk5TunoR75NqKMI2USiXUqOHexdorKsQ/hJa1W7L+4HqzowgfJKVSCYoCReZNrO+UWzvdysz1Mm5FuJ+USiUVFUFurtkpSudn86N3VG8W7FpgdhThY6RUKunmm2HGDDh8+jC/7/zd7DiXdLV6NfO1+RQ5PPy0SngVKZVKiomBQ4cg+1i4xw2CK6EoCje0u4G5m+aaHUX4ECmVKnjoIXj+ORsZx0PMjlKqro27sjZtLTkFOWZHET5CSqUKgoLg+ech+YuBrNM8dxSrrBUk3ElKpYpCQ+G554r4z1MnOHzY7DSXFl07muPZxzmefdzsKMIHSKkYoEvTjjTp9xO//WZ2ktKNSxzHGytNWFJR+BwpFQPYbXaaHBtFnz5mJyldeEg41QOrs/PETrOjCC8npWKQU0fDaNzY7BRlk6kRhDtIqRjA4YDAQLNTlC/YP5h29dqxKnWV2VGEF5NSMYCi6MViBSPaj2D2xtkeNRWm8C5SKgaxWeRP0m6zM6jlIL7f9r3ZUYSXssiPgudTFChyFFFYVGh2lHJd3uJyFuxaILPvC5fwM3qHqqomAs8Xb44E3gNOAXcCQcDM0rY1TTtldB53yc+Hv/b9RWZeJv+M/qfZcco1qvMoPln3CaNjR5sdRXgZw0sFaAn8E/gHsA7ohV4eI4EoYDIQWMq2ZedB7N4djq3vwobqz1qiVDrU78DnGz8nMy+TsIAws+MIL2L42x9N02ZrmpYDxAPJmqZtBjYCzYAITdNm6iTbAAAgAElEQVQ2lbFtWYMGwZJF/mRmmp3Eefcm3stbSW+ZHUN4GZdcU1FVtSXQGNhd/JAdCANyytgOdUUWd3rgAUie19sya+80rt6YQkchaRlpZkcRXkQx+qNFVVWDgTnAKOA2YCn625t+QCawrLRtTdOeLm2/ycnJjpAQz70buMTjb2SS0O9vLotWCQoKMjtOuTLyMvhY+5gJ7SeQk5NjicwlrJYXrJk5KyuLuLg4pfxX6lxxTeVhoC76xdqlwET0s5K7gFzgrTK2yxQTE+OCuMZ6+bFsJjytcGW7IEvkBVidtxrCIeiodTIDpKSkWCovWDNzcnJyhV5veKlomjYNmHbOQ59e8JJbytm2tIb1gmlTMx5IMTuK00bHjubh3x9mdKR8EiSqTsapCALsAXRt3JXVR1abHcXr7Q3ey7yMeXyW/hkLTi+g0FHIvIx5nCg8ceY1OUU5fJv5rWVHPUupCACGtxnOr/t/texfZKuIzI4kISiBZgHNuCz0MlbmrKSI8+cQXpK9hG5B3VAUpy9jeBQpFReoVQuOHrWbHaNCFEXhisZX8OWWL82O4lO6B3enmq3ame3d+bvZk7+HVTmr2Ji70cRklSel4gIDB8Kff1pvQFlc3TiSUpPILfDgtUe83ObczQwMHcigsEFoeRq5Duv9v5BScYH8fDhecIC8wjyzo1TY6NjRfLjmQ7Nj+Cx/xZ9a9lr6r/FHwXpvgaRUXCA/H47lHCf5QMU+ivMErcNbczDzICdzTpodxWecWxydAzuzOGsxP2X+RIRfBAFKgInJKsfwwW+ukpyc7Mhvk292DKc4HPDa2AbQcw4TR/QzO47T7DvtFDYv5Pjp48zfMJ/but1mdqQyleS1Eitm9t/ib/rgN5fpGtzV7AhOqz45hXGvxHGiUVcGDjQ7jXNSbCnEBMdAMGwK3ETDvIZE1og0O1apzuS1ECtmTqZiZ9zy9sdFbDboNWIZ2/8uZLFnLmBYprEJY3ln1TtmxxAWJKXiQp0jOtN1yBrmzYMCi82HFBYQRnTtaNYfXG92FGExUiou1LdpX0ICgrntNpg50+w0FXdrp1v5ZN0nZscQFiOl4kK1gmvRrl47OnWCv/+Go0fNTlQxfjY/+jfvz69//2p2FGEhUipu8tBD8MorZqeouEEtB/Hj9h8tMfeu8AxSKm5SrRrkWW8sHIqiMKL9COZsnGN2FGERUipuUFgIU6fCyJFmJ6mchEYJbDq8iez8bLOjCAuQUnExhwOmT4cbb4T27c1OU3lj4sfwXvJ7ZscQFiCl4mJLl0JGnYW0bWuNkculaV6rOek56RzNstjVZuF2Uioulp0NYdWKSDlqnZngSjMucZzMvi/KJaXiYgMGQGpSF37b+qfZUaqsTkgdagbVZMfxHWZHER5MSsXF7HZ45vEwZr/Slnxr3A9Zprvj7+bd1e+aHUN4MCkVN2jQAOKHrGTqVAcWuSm8VIF+gXRr3I2FuxaaHUV4KCkVN+nSoTaxffdbcgDchYbGDOVb7VsZECcuSUrFTa5rex1DLq9PtWqw3uL36CmKwsiOI5m1fpbZUYQHklJxk2D/YALsAVx5JaxYYXaaqottEMv249vJyM0wO4rwMFIqbpaSAqpqdgpjjEscxxtJb5gdQ3gYKRU3W7kSEhPNTmGMhtUaEmAPYPfJ3WZHER5ESsXNgoMhM9PsFMa5N+Fe3l71ttkxhAeRUnGzYddn854X3UIT7B9MbINY/txr/cF9whhlTnytqmoMMALIBRTADszSNO1vN2TzSh9oT1OvxlOkpECMteY/LtX1ba9n4s8T6dakGzZF/p3ydaX+DVBVdTjwT2CapmlPapr2BPAkcJ2qqv9wV0Bv06haI64asZ933sHyA+FKKIrCje1vlDlXBFD2mcoPmqZlA6iq2hVoVfx4KrDE1cG8Ve+o3vyVupghQ0Ywbx4MHWp2ImN0bdyVeSnzOJ13mtCAULPjCBOVeqZSUijFplzwXI7LEnm5NnXbkHI0hcsug6Qks9MY697Ee3lrldzF7OucXUysDjAWyAYcgAylrCRFObvQW2CgiUFcILJGJEWOIvaf2k/j6o3NjiNM4uxVtRTgVuA2YJTL0viImPAYTued9pprKucalziON5PeNDuGMJGzpeIPPAQ8Bkx3WRofMaLDCEIDQomIgE2bzE5jrLCAMNrUbUNSqpe9txNOK+vTn1tUVe0GoGnaKE3Tbtc0bRQwW1XVIW5L6MXGjIH33oMcL7tCdXOHm5m9YTYObzwVE+Uq9ZqKpmmfqqo6SlXVf6F/4qMADYGNmqbNdVdAb2a3w+TJ8Npr+rpA3sKm2BjeZjhfbvmS69peZ3Yc4WZlXqjVNG0GgKqqtYu3j7sjlC+JjIQML7zRt1dULyb/OpmrWl1FsH+w2XGEGzl1TUXTtONSKK5jt5udwDXGJoyV+4J8kFMfKauqGg10pbiENE2Tj5QNkJSaRGKjRBRFH117zqfNXqFZrWbkFuaSeiqVRtUbmR1HuImzn/58CISiX1fxsr/65pmvzQesuRyqs8Z3Gc/rK183O4ZwI2cHv/0JfKBpWpErw/iaGoE1WLvlFHXqVPe6s5QSYQFhdKjfgb/2/UX3Jt3NjiPcoMwzFVVVF6qqugC4GtitquqC4m1hgJ6RPfnvq4e5916zk7jWTe1vYu6muRQ55N8kX1Depz/93BXEF8U3jOdI7nKCgqLNjuJSJRNlz1w3k1GdZUC2t3Pqmoqqqtef82u5nG8Qf7s/uTmKVw7Xv1B8w3i2H99Oek662VGEi5VbKqqq9gGuUlW1t6qqlwHNXR/Ld9wwOJxffzU7hXtM6DKB11a+ZnYM4WLOnKn0BaKAfkB34FFXBvI1990Yw+LFZqdwj/ph9akVVIutR7eaHUW4ULmf/mia9riqqpoMzXcNm817B79dypj4MTz020O8fMXL500DIbyHs+NUklVVfVFV1edVVZW3PwarVg3SfeRSQ4A9gH+0+Ac/bP/B7CjCRZwtlSeB2cDXgIxkMtiAAfDTT2ancJ+BLQfyx84/yCv04lF/PszZUknTNG0tkAzsd2Een9S5MyxfDoU+tN75mPgxvLfai9YqEWc4Wyr9VFVNBlYBCTIAzlhFjkJqd/+G133oHLB1eGuOZR/jyOkjZkcRBnO2VBKAm4CrgK6apl3muki+x26zk193NTVqwKJFZqdxnwldJsh9QV7I2VL5BHgCaAt85LI0PizAHsBNt+TyzTfesx5QeWoF16Jx9cZsOLTB7CjCQM6Wykngb03TfnZlGF+W2CiR1Wmr6NQJ/vah9R/viL2Dj9Z8JFNPehFnS+UYMEBV1Q+Bwy7M47O6N+nOn3v/pFs3/aKtr/Cz+fGvVv/iW+1bs6MIg5S3lnJH4FkgADgFRAIvuiGXz6keWJ2MvAxatYJZPjYF1uUtLuf+X+5nYPRAAv28bDGkC3y+43NWrF5BWmYaXRt15XT+aU7lniI0IJTPh33OkdNHGPnNSGoE1uDxvo8TU9d6C26Xd6byBHCrpmn9NU0bgH6x9knXx/JNj/R6xGvnVSnPPfH38O7qd82O4XI3Rt/If3r9h6taXcXkHpMZ0X4E39/0PaH+oWw5soXpi6bzTP9n+PSaT3lp+Utmx62U8obpFwJBqqpGnvOYs2+ZRAWF+IcA4Ofs1FlepGWdlqRvSudQ5iHqh9U3O45btKrTilZ1WpFfmM+uk7toVacV+0/tJ7ZBLABBfkEmJ6yc8v76pgOPX/DYKRdlEcV89WxlfJfxvPDnCzzd/2mzo7jV1IVTmdhlIgH2gDP/sACcyDlhYqrKK2+SJplRxwRFRd45EXZ5agbVpFmtZqw7uI5OEZ3MjuMWM9bOwN/mz9WtrwagY/2O/LXvL9rVa0duQa7J6SpH3sp4oIYN4cABs1OYY1SnUcxYO8MnPmLeeWInE36ewMHMg9z13V1sOryJexPv5Z3V7zBy3kge7vmw2RErRbHK/7zk5GTHdxnfmR3Dab1CerE0a2mFv6/IUcSu9Y0hP5gWCe4dsFLZzEbbdWIX2QXZtKnbpszXeUreirBi5quqXUVcXJzz580Oh8MSX6tXr3ZYyZYtWyr1fUdOH3E89fPbjldeMTiQEyqb2RXu//l+R1ZeVpmv8aS8zrJi5uKfPad/VuXtj4cJDwknx/8AJ6x5jc4w4xLH8daqt8yOISpBSsUDKYpCQYE13pa6SrNazSgsKmTPyT1mRxEVJKXigRopsRSF+uiV2nPc1+U+3kx60+wYooKkVDxQ9ZO9yQlfYXYM04X4h9ClcRcW7JLpe6xESsUDpe2sTXRra45RMNqwmGF8u/VbCooKzI4inCSl4oGOH4d7e99kdgyPoCgKd8bdyfvJ75sdRThJSsXDFMlywxdpV68dhzIPydSTFiGl4mE2bYL27c1O4Xkmdp3IqyteNTuGcIKUioepX9+3Zn5zVq3gWrSo3YLVB1abHUWUQ0rFw9SvD1FRMHu22Uk8z60db2XW+lkUOeQ9oieTUvFAI0boqxbeN/kEOw6mmR3HY9htdm5sdyOzN0jjejIpFQ81eDDccOsp7ptyhFmzIE8W8wOgW5NubD6ymVO5Mq2Pp5JS8WDd20YSe8sXtGsHjz0GTz0Fycm+s4RHaSZ0mcBrK14zO4YohQ9OXGgdiqKgKAqxsRAbC5mZ8PPP8PXX+gROzZtDnz76f31pQqcG1RpQI6gGO0/tJAbrTQzt7aRUPFxkjUj2nNxDVM0owsJg+HD9y+GAXbv0FQ0//lh/bdOm0LcvREd7f8ncHX83o78YzZWJV6J4+2/WYqRUPFyfqD4s2bOEW2rect7jJWcqzZuffWz3bli8GGbO1EsnLg6uugr8/d2b2R0C7AH0atCL77Z9x2B1sNlxxDmkVDxcqzqtyMzLdOq1TZvqX6CXysqVMGkSTJkCjRu7LKJpekb05MPdH/KPFv+w7Mzz3kgu1Ho4RVGIaxhXie+Drl2hXTtITXVBMA9xb+K9vJUkkzl5EikVL7V9O0yeDI0aQZcuZqdxnea1mpNflM++9H1mRxHFDH/7o6pqDeA19KVS3wUeRl8r6E4gCJhZ2ramaTL4oAqysvRPh1au1N8GPfkkBPnAu4L7Eu/j8cWP8/zlz5sdReCaayoJwFSgBrABaI9eHiOBKGAyEFjKtkzzVUF//w2//AL790NICFx+OTz7rPd/+nOu0IBQ4hvGs3j3Yvo07WN2HJ9neKlomva7qqpRwBRgnqZpm1VVDUBfh7mepmmbytgW5Th9Wv8YecUKKCyEFi300bfeeCG2Iq5tcy0Tf55Ij8ge+Nnk8wczuepP/xiwGbi9eNsOhAE5ZWyHlrfTlJQUY1O6UE5OjmF5i4rg+dm7KTzYjjrVAklMzOK667LPrLmckQFGHMrIzO5wYd5+Nfvx5A9PckP0DSamKpvV/owrwxXXVIYBi4DPgcdUVY1Ff/uzF8hUVTUO/e3OhdvlXmmLibHO6MmUlBTD8r77LvSI9aN2m5UMbzPckH1eipGZ3eHCvDHEsHbhWupG1SU8JNzEZKWz2p8xQHJycoVe74pPf7YCs4AXgXhgEvp1k/eBGcDEMrbFJcTHw6GtLVh/cIPZUTzehK4TZDInk7nimspmYNA5D91ywUvK2xYXiI+HggIb/36nL45+vnURtqJqB9emac2mrElbQ2yDWLPj+CQZp2IRXbtCTOxxZs3NMDuKxxvVaRSfrPvEJxZ590RSKhYyengUX/96iMJCs5N4NrvNzvVtr+fzTZ+bHcUnSalYSFzDOP59V1PmzjU7iefrEdmDDYc2kJErZ3buJqViITbFRo9ufqxdC7my1li5xncZz+srXzc7hs+RUrGg/v1htUwqX66G1RoSFhDGtmPbzI7iU6RULCgqSh+WL8p3T8I9vLPqHblo60ZSKhZUqxacOGF2CmsIsAdwRfQVfL/te7Oj+AwpFQs6mXNCZtevgH9G/5NFuxeRU5BT/otFlUmpWIzDAXc9onHV0NNmR7EUmczJfaRULOb552HEDf7syP3L7CiW0rxWcwqKCtibvtfsKF5PSsVCPv5YX7z91ivbsDJ1pdlxLOe+Lvfxxso3zI7h9aRULOK778BmgyuvhGD/YLk+UAkh/iF0bdyVP3b+YXYUryalYgErVsDWrXDbbWcfC/ILIjs/27RMVjU0ZijfbfuO/MJ8s6N4LSkVC/j8c3jwwfMfG9B8AKkZXjxNvosoisKYuDG8l/ye2VG8lpSKBfj767O7natr465E1442J5DFxdSN4Xj2cQ5mHjQ7ileSUrGAhx+GRx6R+32MNLHrRFnk3UWkVCygTh3497/1cpFpD4xRPbA6beu15a998tG80aRULCIyEkaPhunT9QFwoupGtB/B3E1zKSySpjaSlIqFtGkDgwbBCy+YncQ7KIrCqE6jmLFuhtlRvIqUisV066YPgHuv+MOLxbsXmxvI4jo36Myek3s4lnXM7CheQ0rFggYOhIAAWLYMFu5eSEFRgdmRLG1i14m8tlIu2hpFSsWibr4ZFi6E2AaxrElbY3YcS6sTUocm1ZuwNm2t2VG8gpSKRfn7658E9YzsydI9S82OY3m3d76dGetmyGROBpBFZy3M4dDXuTmRIzM2VZXdZueGdjcwZ+McRnQY4bLjfPvtZh5/fBGHD2eQmNiMyMg67Nx5BEWB5567lmPHTjN+/GyqVw/iwQcH0rJlfZdlcRU5U7Gwkn9UbYpNPhY1QPcm3dl4eCOnck+57BhXX92WCRMu5/LL2zJt2mBq1w7lxRevJzQ0kK1bD/Liiz/x8MODeP31m3n33YUuy+FKUioWVr067NsHN3e4mbxCmQrOCO4caRsSEsgNN3Thiy+SSEtLR1UjSEtLp0OHJgQG+hEY6O+WHEaTUrGwsWP1MSs5qa0I9g82O45XiAiLoGZQTbYe3eq2Y7Zv35js7DwOHTpFcPDZIklPz3JbBiNJqVhYSAi88gosWQIvvgj5cje/Ie6Ov5t3V7/r8ou2+/ef4Lvv1tGmTUM6dmzCvn3HadOmEatW7SIjI4e8PGsOFZALtRZnt8O4cbBtGzzwAIwaBZ07m53K2vzt/gxWB/N1ytcMazPMZcdp2LAGzz+/mZ9/3khISAAJCU1p3boBjz76FZmZuUya9A+XHduVFKt8hJacnOxYuPCI2TGcFhdXk+Tkk249ZlERLF8eicOh0L37HmwVPA81I3NVuDrv93kfMcD/RoKUEMP2abU/Y4B+/eoSFxenOP0NDofDEl+rV692WMmWLVtMO/aaNQ7HuHEOR2pqxb7PzMyV4eq8B04dcDyx6AlD92m1P2OHw+Eo/tlz+mdVrql4kZJJnTt3hmeegQ8+gA8/1M9gRMU1qNaAOiF12HR4k9lRLEVKxYscOn3ozMXFsDB9moSOHWHCBDhwwORwFnVX3F18kPyBjLStACkVL9KqTquLFiNPSNDXCnrmGTljqQw/mx83tr+RTzd8anYUy5BS8SJ9ovqweM/FUyEEB5+d4Ek+dq64ro27su3YNo5nHzc7iiVIqXiRqJpRpa7A17GjfmfzI4/ArFlQYM0hEKaZ1HUSryx/xewYliCl4mVsSun/S1VVH4Hbrh089BB89ZUbg1lcnZA6tKjdgqTUJLOjeDwpFS8zvc/0cl8TGwsvvQR5efD7724I5SVGdhzJ7A2z5ebNckipeBm7ze70a8PDqfAAOV9mU2zcEXsHH6750OwoHk3+SvmwxYuhd2+zU1hLh/odOHT6kCxEVgYpFR9mt4Of3P1VYRO7TuTVFa+aHcNjSan4MPkEqHKqB1YntkEsi3YvMjuKR5JS8UJFjiKOnC775suvvoLERDcF8kLXtrmWb7Z+I5NjXYKUihdSUMpccmLTJn2qhKuvdmMoL6MoCmMTxvL2qrfNjuJxpFS8kKKUfZf6Z5/pazOLqmlVpxW5BbnsObnH7CgeRS7TeakGYQ1Iy0ijQbUGFz2XkKAv9l67NgQF6V/BwXDsWDW2bdO3GzeGtm1NCG4x47uM59EFj/LSFS+ZHcVjSKl4qd5RvVmyZwnXt7v+oueGDYNrrtEHv+Xk6F/Z2bB5cy6NGunbb78Nr76qry8kShfsH0z/5v35ftv3/KvVv8yO4xHk7Y+XaluvbZnzgNhs+hlJzZoQEQHNmkGLFnl07qyv1zxwICQnuzGwhV3Z8koW7lpIVr41J6o2mpSKl7IpNtrVa1fp7+/bF2bOhLQ04zJ5swldJ/D6ytfNjuERpFS82KXe+jgrLEyfh2XWLPjPf2CNLNdcpsgakQT7BaMd1cyOYjopFVGqatX0u5mnTdM/gv73v/XpKTMyzE7mme5JuId3Vr/j87PESamIcgUGwg036GcuvXvDyy/Da6+dXXZV6ALsAVzT+hr+t/l/ZkcxlZSKqJBWrfQZ5OLjYcYMs9N4nj5N+7AmbY1L12P2dFIqolJ69IDt281O4ZkmdfPtWeKkVLzcsaxjrElzzVXWOnXgpLXWxXKLiLAI6ofVZ93BdWZHMYWUiperHlid+dp8l+w7NFQu2pbmztg7+XjtxxQ5fG8JAykVL+dv96egyPg5DvLz9bc/TZoYvmuvYLfZGdlxJDPW+t6FJykVHxAWEEZGrrGnFB98oC/7IUoX3zCePel7OHz6sNlR3EpKxQf0aNKDP/f9adj+Cgth925o08awXXqtB7o94HMXbaVUfEBio0RDl5aw2/UbELOzDdul16oRVIOOER19apY4KRUfEOgXyPA2ww3d58MP6yNsc3IM3a1Xur7t9cxLmeczs8RJqfiINnWNfa/SsKE+hP+RR2RkbXkURWFc4jjeTHrT7ChuIaUiKq1xY+jZE5YvNzuJ52tZpyUFRQXsy9xndhSXk1IRVXLVVfDDD2ansIbxXcbz2fbPvP6GQykVUSX+/lC9Ohw4YHYSzxfkF0SfBn34OuVrs6O4lJSKqLKxY+G998xOYQ3dI7qzMnWl4eOGPImUig9ZuX8lf+41brxKiYMHISTE8N16rfu73c/Ly182O4bLSKn4kLb12rJ071JD97l6tT5x0wMPGLpbr1Zyw+HatLVmR3EJKRUfEhYQRmZepmH7++47+OMPePZZWZO5orz5hkMpFR/jb/M3ZBDWJ5/AiRP6WJVy1i4Tl2C32bmt0218vPZjs6MYTkrFx8Q3jCf5QNXW3ti+HdLTYeRIg0L5qLiGcexL38ehzENmRzGUlIqP6RnZkxX7V1RpH2vW6DO/iarzxou2Uio+pkZQDSZ2nVilfQwZAnPmyEhaI9QIqkFcwzgW7FpgdhTDSKn4oPIWcC9PYCC89BJs2QJPPSU3FVbVtW2uZb42n9yCXLOjGEJKRVSKosAdd8DNN8PkybBjh9mJrKvkhsM3kt4wO4ohpFRElTRtqq8D9LH3fYjhVtG1o3E4HOw4bv12llIRVebvr3+JqhnfZTxvJr1p+RsOpVR8VFpGGoVFhYbtT1H0aSZF5QX6BTJYHcz/bfk/s6NUiZSKj1q2dxkbD280bH+DBsHs2Ybtzmdd1uwyVh1YZekVDqVUfFSvqF4s2bPEsP0lJOiLuB/yrnFcpnig2wO89NdLZseoNCkVHxURFmH40hEPPaQv4m7xSwKmqx9Wn8bVG7tsZUlXk1LxcUZeFKxWDYYPhyeegALj1y/zKXfE3sGMtTMMve7lLi67t1RV1S+BF4DHgVPAnUAQMLO0bU3TrPtG0oLa1G1DytEUQyfF7tYNatSASZPgnntkbaDKsik2bu98Ox+u+ZAx8WPMjlMhLjlTUVV1AJAAPAhMBp4HRhb/uqxt4Ub9mvYjPSfd8P22aQOvvAK//66/Hcr1joGibte5QWfSMtM4mHnQ7CgVYnipqKrqB9yGfgaSq2naJmAj0AyIKGdbuFGDag3o1qSbS/bt5wfjx+tvh/79b/jrL5ccxus90O0By91w6Iq3P+OA14GBQIPix+xAGJBTxnZoeTtOSUkxNKgr5eTkWCovuC7zmDHw+ee1WLGigIEDjZub1Vf+jBs6GvLRoo/oXr+7i1IZyxWl0gmIBLoB8aqqxqJfO9kLZKqqGgcEXmK73AVRYmJiXBDXNVJSUiyVF1yb+ckn4fPP9WkTRowwZp++8mfcunVrJv0yiREtRxDkF+SiZKVLTq7Y/DuGv/3RNO02TdPuB34CYoFJ6NdN3gdmABPL2BZe7MYb9cmdVq40O4m1KIrC+C7jeX3l62ZHcYrLPv3RNO2J4l/ecsFT5W0LL3b33fok2V26mJ3EWprXao5dsbP92HZa1mlpdpwyyTgVwU/bf3LbTWw2G/TqBUuMG8zrM0rWY/b0Gw6lVASpGansPrnbbccbMgS+/lpG3lZUoF8gQ2OGMnfTXLOjlElKRdAnqo+h9wGVx2bTB8mt9c5lb1yqT9M+bDi0gePZx82OUiopFUF07Wi3Tw60axc0auTWQ3qNB7s/6NE3HEqpiCrPWVtRhw/r9wbVr+/Ww3qNOiF1aB3e2iVL2BpBSkUA0Kh6Iw5muWc4+IoVMGCAWw7ltW7ucDNzN801ZGE4o0mpCACGtxlOsF+wW47VqxcsNXZJZ59TMln2W0lvmR3lIlIqAoDwkHBqBNRwy7Hsdn0QnKgaNVwlvyifnSd2mh3lPFIqwq3S0+Hhh+H++81O4h3GdxnPayte86ixK1Iqwm3S0+GRR/T7gGrXNjuNdwjyC2JI6yH8b/P/zI5yhpSKcJtffoG77pJCMVq/Zv1Ye3AtJ3NOmh0FkFIRbhQcDPn5ZqfwTg90e4AX/3rR7BiAlIq4wHPLnnPJfouKYNkyfUVDYby6oXWJrh3N8n3LzY4ipSLOF2AP4GjWUUP3mZamz1l7/fUQHm7orsU5bu14K3M2ziG/0NzTQSkVcZ5eUb1YtneZYfv7+mt4+234738hNtaw3YpLUBSFsQljeWuVuWNXpFTEeTpFdGJtmjF3+r35pn7z4JNPQmi5k4UKI8TUjSGnIHYJLfsAAAtHSURBVIddJ3aZlkFKRZzHz+ZHoaPqa83MnQsNG+rTHAj3mth1Iq+tNG/sipSKuEh4SDin805X+vuTkyE1FYYONTCUcFqQXxD/avUv0xZ6l1IRF5nYdSKhAZV/v/K//8GECQYGEhU2oPkAVh1Y5ZJ1ncojpSIMd8cd+qc969ebncS3Pdj9QV5a7v55V6RUhOFatYKXX9ZL5d//hvfeg4PWWmTPK9QLrUfTmk1Zud+9yxe4bDZ94dv8/WFk8UK2u3bBF1/oxRIWBj16QNeuEOT+JWx8zm2dbmPCTxOIaxiHn809P+5SKsLlmjU7e43l9GlYvhxefFFfY9luh9atIT4eWrQAN09C5/Vsio17Eu7hnVXvcF+X+9xyTCkVUardJ3fTtGZTQ/cZGqrP+lYy81thIWzdCn/+CbNm6TPsBwToA+W6doU6dQw9vE9qU7cN32z9hr3pe4msEeny40mpiFJ9tOYjnrzsSZcew26Htm31rxK5ufryqLNmwfHjetG0agWDBknJVNbErhN5dMGjvHyF6xd7l1IRpQryCyIrP4sQ/xC3HjcwUF/Co1s3fdvh0M9mZsyAI0f0sxgZA1MxIf4hDGg+gPnafAarg116LPn0R5Sqa+OurNi/wuwYKArExMCDD8Jzz+nXaB59FD77rBbZ2Wans44rW17Jkj1LqjSw0RlSKqJU3Zp084hb6S+UmKiXS8+emTzxBLzxBmRlmZ3KGiZ2ncgrK15x6TGkVESpQvxDyC7w3FOBpk3zeeYZGDgQnnoKPvxQllItT+PqjakVVItNhze57BhSKqJMiY0SzY5QruhofWqF9u31kbwyU3/ZxsSP4f3k9ylyFLlk/1IqokyuvqhnpC5dYOpU+M9/YNs2s9N4Lj+bHyPaj2DW+lku2b+UivAqderAq6/Cd9/pS4H8+KNcb7mULo27sOP4DsNn+QP5SFl4IT8/eOABfb3mlSv1+5Cys6FePbj8cv2TJBm5q0+W/cJfL/Df/v81dL9SKsJr+fnp9xn16KFvHzoEv/8Os2frpdKhA/Tv77sD6moF16Jt3bYs3bOUXlG9DNuvlIrwGfXrw4gR+q8dDtiwAT7+GI4ehfHjoVEjc/OZ4ab2NzHh5wl0adyFAHuAIfuUayqiXEWOIr7f9r3ZMQylKNCxI0yeDI89Bk8/7ZtrEimKwr0J9/L2qrcN26eUiiiXTbF55CA4owQH62+DNm82O4k51HCVnIIcdp/cbcj+pFSEU+w2O4VFVZ8Q21Pt3Qs1a5qdwjwTukzgjZVvGLIvKRXhlE4RnVh3cJ3ZMQx38CBMmwYREb69emKwfzD9mvXjh20/VHlfUirCKb2jerN071KzYxjm2DF9aP+MGXD//XDjjWYnMt+/Wv2LP3b9QXZ+1W7NkFIRTgkPCedY1jGzYxjC4YApU2DMGH2AnC+/7bnQ+C7jeSOpam+D5CNl4bQ7Yu8wO4IhFAXefVefIEqcr2nNpoT4h3Aq9xTVA6tXah9SKsJpRk8taSYplNKNSxxXpe+Xtz9CCENJqQghDCWlIoQwlJSKEMJQUipCCENJqQghDCWlIoQwlJSKEMJQUipCCENJqQghDCWlIoQwlJSKEMJQUipCCENJqQghDCWlIoQwlJSKEMJQUipCCENJqQghDCWlIoQwlJSKEMJQUipCCENJqQghDCWlIoQwlJSKEMJQUipCCENJqQghDCWlIoQwlJSKEMJQUipCCENJqQghDCWlIoQwlJSKEMJQUipCCENJqQghDCWlIoQwlJSKEMJQUipCCENJqQghDCWlIoQwlJ/RO1RVtQswHdgNfAI8DpwC7gSCgJmlbWuadsroPEII9zK8VIBOwGRN0zarqvoCMBkIBEYCUeVsv+mCPEIIN3LF25/OwCRVVb8BOmmatgnYCDQDIsrZFkJYnCvOVB7UNC2z+G3QsuLH7EAYkFPGdmh5O05OTjY4qmtZLS9YL7PV8oI1M1eEK0rlJuB9oDYwTVXVOPS3N3uBzDK295W107i4OMUFWYUQBnNFqeQUv/UpACYB/0U/K7kLyAXeKmNbCGFxisPhMDuDEMKLyDgVIYShpFSEEIZyxTUVQ6mqWg8LDJBTVTUReL54cyTwHh6euYSqql8CL3DOQEVPzqyqajvgZuBlPPzvRvGnoBPQP928Dw//e6Gq6i1ATeALKjlQ1eNLBX1wnBUGyLUE/gn8A1gH9EL/H+HJmVFVdQCQADyINf6cAR4D7gAexfMzjwDuBToAa4HeeOjfC1VVOwB9gA1c/HPn9EBVK7z9scQAOU3TZmualgP8f3tnF2JVFYbh5zBKTGJlECVFRWZfP8qQDmSl1Fx0EYQ3GTbDEDMSWGHmRIgaNpWWRDf9wNCfTPQnaYR1UUrJyIiRQtTEjPLSDxbaRUrpKIaVc7r41qHDwTE9s6c5m/meq7P2OWvtd2/2+vZae3/vOo3AV5IGqHHNZjYBaMPvQCfycJ7NbAFwPdANXJUDzVvxDrgU2FLL14Wkb4G3gAIjSFTNQ1A5qwS5scTMpgOX4b4nqH3NS4CXgCFgatpW65pbgDvxIfmUtK2WNc/BpxIfADenbbWst5QPVnWiah6mP/1nmiA3lphZPf5MpR1oy4Nm3Kd1OX6xN5rZLHxoXsuafwcO4Od2dg40zwDWAOcBXTnQW6K/2kTVPASVbvKRILcSuAgPLDuAZdS4ZkltAGb2BK5zOTWuGXgFeBuYDNyFJ1jWsuY3ge34aHA+ta8XoIj3uy6qSFSN5LcgCDIlD89UgiDIERFUgiDIlAgqQRBkSgSVIAgyJYJKEASZEkElyD1mVmdmzVXWvdHMZmataTwTr5Rzipl1ArcCR4ABSU+OrSLHzK4AtgFf4BmviyQdHOV9rgA2Aqtxw9tTwCHgSmA38DlukntE0g8VdQvAC/hi7X+Ops7xQoxU8s1aSfcAs8xs8liLKaNX0n1ADzB3NHeUMpnrJP2YNk0HGoA/UnmPpFb8b2MeqqwvqYgned07mjrHE3nIqA3+m3OA181MwDRJrWlI34H7NN4FPgHew+/gxyQtN7OH8X8/mALcj3uBeiT1mlmPpCYz6wE+BaZK6jCzJtx5eyHwnKRdp9Azz8w24M7cRjObh9v/h/BRQTvwIrABWJvqFIDbgUuAR4GTuBHvEPCRpM1mth73oNRLWpTqXYe7f0sU036/SW2WvCznA4NmdnXa9zHgQ0nvS5KZLTzz0x2cjhip5JuVZvYO8AZ+Z96f7soAjwOHcd/GLbjf5FI8QLyaftOAd8g1eCcrp3xe/LWkjvS5EzgI7Odfg1wlOyQ14523BZ+WNONroHQCX6byPtx/tAtYBfwG/IovxQBwDfCgpM2pbHjK+/MV+ytfFP1noAn4LpWvNbNuYCEe0AbT8e8FFgxzvMEIiKCSb9ZJapW0KZW3ln1XxBcxWo3b1cFHIgXcOwPu5dgNrABuwkcSE81sEnBxqSFJn5W1exJfz+RZ4PthdJU6+WHg3NRuAb/eingQaQY2AbdJ2gcclbQqaT6Q6u9My0mUlmlYltpcnzQCDACzy/bdl7QPpfJeSe2SFks6grubG3DX8ITU9ozTHEtwlsT0J79U3lkry+uAl9Pn19L3jwFHgf60/Q5gGlCPjzx68cAzJ/3uVO0+g09b6oCnh9E2N42gLsDXa+nDp17FVGcPcBw3Xt6d6nSZ2UZ8Kre0skFJf5vZA8BfwC+kZyaSTpjZoJndkNrvwxcZmjiMtp/wxbSWAJNSsFqMTxWDDIi3P8GIMLMW3LlawDv1zrQQ0f+tY76kj6uoNxM4XvlWKKieCCpBEGRKPFMJgiBTIqgEQZApEVSCIMiUCCpBEGRKBJUgCDIlgkoQBJnyD8puQ+KpKi38AAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x1295b4e0>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "fig_pres, ax_pres = plt.subplots()\n",
    "ax_pres.invert_yaxis()\n",
    "# plot hydrostatic\n",
    "well_cug1.hydro_log().plot(ax_pres, linestyle='--', color='green', label='Hydrostatic')\n",
    "# plot OBP\n",
    "well_cug1.get_log(\"Overburden_Pressure\").plot(ax_pres, color='green', label='Lithostatic')\n",
    "# plot pressure\n",
    "pres_log.plot(ax_pres, label='Bowers', color='blue')\n",
    "# plot horizon\n",
    "well_cug1.plot_horizons(ax_pres)\n",
    "\n",
    "\n",
    "# set fig style\n",
    "ax_pres.set(ylim=(5000,0), aspect=(100/5000)*2)\n",
    "ax_pres.legend()\n",
    "fig_pres.set_figheight(8)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.15"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
